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ABSTRACT 

Riemer, Terry Edmund, Ph.D., Purdue University, December, 
1974. Optimum Constrained Image Restoration Filters. Major 
Professor: C.D. McGillem 

An optimum image restoration filter is developed in 
Hilbert space by minimizing the radius of gyration of the 
overall or composite system point-spread function subject to 
constraints on the radius of gyration of the restoration 
filter point-spread function, the total noise power in the 
restored image, and the shape of the composite system fre- 
quency spectrum. 

By satisfying a constraint on the radius of gyration of 
the restoration filter point-spread function, truncation 
errors arising from attempted restoration of incomplete 
segments of blurred data are almost completely suppressed. 

An iterative technique is introduced which alters the 
shape of the optimum composite system point-spread function, 
producing a suboptimal restoration filter which suppresses 
undesirable secondary oscillations which may otherwise 
appear in the composite system point-spread function and 
introduce w ghosts” in the restored data. 

An extensive study of the restoration filter performance 
as a function of its parameter variations is made. Numerous 
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examples are provided to demonstrate the fundamental prop- 
erties of the restoration filter. 

Finally this technique is applied to multispectral 
scanner data obtained from the Earth Resources Technology 
Satellite to provide resolution enhancement. An experimental 
approach to the problems involving estimation of the effective 
scanner aperture and "matching" the ERTS data to available 
restoration functions Is presented. 


CHAPTER 1 


INTRODUCTION 


1.1 General Discussion 

The general problem of image processing has been the 
subject of intense investigation within the last decade. This 
interest is the result of the need for the highest possible 
image quality in the increasing application of imagery to 
the solution of various related problems in many fields of 
science and engineering. These applications range from x-rays 
in medicine to satellite-based multispectral scanner data 
collection systems for monitoring earth resources and to a 
large degree are made feasible by recent improvements in 
digital computer hardware. 

Since no image collecting or imaging system will produce 
a perfect replica of the original image, some further pro- 
cessing is usually required. The area of image restoration 
deals primarily with the problem of processing the output of 
an imaging system in such a way that the significant parameters 
or features of the original image are, in some sense, enhanced 
or restored. This processing may be linear or non-linear, 
shift-variant or shift-invariant depending upon the type of 
degradation produced by the imaging system. 

The primary purpose of this research is the development 
of a technique to reduce the effective aperture radius of the 
multispectral scanner system used for remote sensing of 
earth resources. Because of the finite point-spread function 
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of the scanner and the response limitations of the image 
sensors and signal conditioning electronics, a two-dimensional 
spatial smearing or blurring of the original image is pro- 
duced. This type of imaging degradation essentially maps 
many points from the original image into a single resolution 
element. In other words, a single resolution element of the 
imaging system output represents a two-dimensional weighted 

l 

sum of many points adjacent to the correspondingly located 
sample of the original image. Thus, depending upon the 
density and shape of the effective aperture and the spatial 
and multispectral characteristics of the original image, 
this type of imaging degradation may not only introduce 
significant error in the visual presentation of this data, 
but also introduce serious errors when classification of 
this data is attempted. 

1.2 Previous Investigations 

Historically much of the initial work directly related 
to linear shift or time-invariant restoration filters was 
performed in the areas of seismographic data analysis and 
radio astronomy. Basically this work considered the prob- 
lems associated with defining the restoration filter trans- 
fer function as the inverse of the Fourier transform of a 
linear shift or time-invariant blurring process, which is 
often referred to as a classical "inverse" filter. The 
principle problem associated with such a restoration filter 
occurs when the Fourier transform of the blurring process 
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becomes very small or perhaps equals zero. Since the 
spectrum of the restoration filter becomes unbounded at such 
points, some means of assigning a finite value to the 
spectrum at these points is necessary. 

Typical of this early work, Harris [9] attempted to 
eliminate this difficulty by arbitrarily modifying the 
spectrum of the inverse filter to provide a restoration 
filter having a transfer function that was continuous and 
finite in those regions where the blurring function spectrum 
is zero. Bracewell offers an alternative to the problem by 
employing finite differences [4] and a technique of suc- 
cessive substitutions [5] to define the restoration filter 
in the temporal or spatial domain. In a more recent study, 
Kinzly [13] circumvents the problem of infinite gain in the 
restoration filter by choosing a diffraction-limited system 
as the correction goal rather than attempting an ideal 
correction. More sophisticated spectral studies by Bercherer 
[2] and Horner [10] are based upon minimizing the mean-square 
error between the original and restored image data, leading 
to a two-dimensional Wiener filter. Although this type of 
restoration filter does eliminate the previous problem of 
spectral infinities, "Apriori" knowledge of the image power 
spectral density is required. Using a novel approach, 
Robinson [18] defines the restoration filter in tne temporal 
or spatial domain in terms of a series derived from the z- 
transform of the blurring function. Rice [17] provides an 
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extension to Robinson’s work by defining a restoration filter 
based upon minimizing the mean-square error between the 
original and restored images. In both of these procedures 
knowledge of signal characteristics is required. Extensive 
literature surveys on restoration filter design and related 
areas of image processing are provided by Huang, et.al. [11] 
and Sondhi [20] . The previous references represent the 
earlier studies most closely related to this investigation. 

Although several of the previous studies did attempt to 
optimize the restoration filter, resulting in a Wiener filter. 
Smith [19] appears to be the first to consider specifying an 
optimum restoration filter based upon minimizing the radius 
of gyration of the corrected impulse response or point- 
spread function subject to a noise constraint. Because no 
constraint is placed upon the radius of gyration of the 
restoration filter impulse or point-spread function, serious 
numerical problems arise. Stuller j.1,21] uses basically 
the same constraints as Smith, but formulates the problem in 
the spatial domain. By defining the restoration filter using 
matrix notation, a restoration filter containing an arbitrary 
number of points may be determined. This appears to be the 
first attempt to control the duration of the restoration 
filter impulse or point-spread function. Although the 
duration of the restoration filter impulse response or point- 
spread function can be controlled, no control is possible 
over the shape of this function. 
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The previously cited investigations have in common one 
or more of the following deficiencies: no general control 

over the shape or duration of the restoration filter impulse 
response or point-spread function, no general control over 
the shape or duration of the corrected impulse response or 
point-spread function, no control of the spurious responses 
resulting from secondary lobes in the restoration filter 
impulse response or point-spread function, and no way to 
make tradeoffs among the various conflicting performance 
requirements of the filter. The significance of these 
controls will be examined in the next chapter. 

1.3 Outline of Investigation 

An optimum shift- invariant restoration filter is formu- 
lated in Chapter 2 based upon minimizing the radius of 
gyration of the overall or composite system point-spread 
function subject primarily to constraints on the radius of 
gyration of the restoration filter point-spread function and 
on the total noise in the restored image. In addition, an 
iterative technique is developed which suppresses secondary 
oscillations in the composite system point- spread function. 

For numerical solution convenience, the problem is transformed 
into the spatial frequency domain and a system of linear 
differential equations which specify the spectrum of the 
restoration filter is developed. Finally a comparison is 
made between the optimum restoration filter and the classical 
"inverse” filter. 
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Chapter 3 discusses the application of well known digital 
computer algorithms to the solution of the previously devel- 
oped system of linear differential equations specifying the 
spectrum of the restoration filter. 

Chapter 4 presents an investigation of the fundamental 
properties of the restoration filter. A study of the filter 
performance as a function of its parameter variations is made. 
In addition, examples are presented to demonstrate the capa- 
bility of the restoration filter to suppress truncation 
error. These examples also provide a qualitative measure of 
the restoration capabilities of the filter. 

Chapter 5 discusses the application of this restoration 
technique to ERTS multispectral scanner data. The problems 
of estimating the shape of the scanner aperture and "matching" 
the ERTS data to previously computed restoration functions are 
considered. Finally, examples are presented to demonstrate 
the effectiveness of this technique in the ERTS data. 

Chapter 6 presents conclusions made in this study and 
suggests related work for further investigation. 
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CHAPTER 2 

OPTIMUM IMAGE RESTORATION FILTER 
2.1 Introduction 

In this chapter the fundamental equations underlying 
the general image restoration filter are defined. From 
these equations, a system of equations which specifies the 
general two-dimensional restoration filter is formulated 
in the two-dimensional spatial frequency domain. To 
simplify the analysis of the general problem, two specific 
types of aperture functions are considered, separable and 
radially symmetric. A system of equations which specifies 
the spectrum of the restoration filter for each aperture 
class is developed. In order to provide additional control 
over the shape of the point-spread function of the 
composite imaging system, an iterative technique is 
defined. The resulting systems of equations which specify 
the frequency spectrum of the restoration filter for both 
separable and radially symmetric aperture classes are 
formulated in terms of the iterative technique. And 
finally, it is shown that when all defined constraints are 
removed, the resulting restoration filter is the classical 
"inverse" filter. 
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2.2 Formulation of Image Restoration Problem 

The basic components of the composite imaging system 
are defined in Figure 2-1. The fundamental design 
criterion for determining the optimum image restoration 
filter is to select the restoration filter point-spread 
function so that the point-spread function of the composite 
imaging system is as compact as possible about the spatial 
origin. Since the classical "inverse" restoration filter, 
defined by the inverse of the modulation transfer function 
of the imaging system, would produce an impulse for the 
composite imaging system point-spread function, 
representing theoretically the best possible correction, 
the basis for the fundamental criterion is evident. 

One convenient measure of compactness is the square 
of the radius of gyration, r, defined as 

00 _ 2 a _ _ 

/ | v | g (v) dv 

2 —oo 

r - (2-1) 

“ 2 „ _ 

J g (v) dv 

— CO 

where g(v) i3 the point-spread function for the composite 
system and v is a two-dimensional spatial vector. g(v) may 
be expressed as a two-dimensional convolution in terms of 
the imaging system point-spread function, h^(v), and the 
desired image restoration point-spread function, h r (v) as 

g(v) = h b (v)**h r (v) . (2-2) 

An optimum restoration filter based upon minimizing 
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Eg. 2-1 would only assure that the composite imaging point- 

spread function be as compact as possible but would not 

provide control over oscillations or secondary lobes that 

might occur and which would lead to ghosts in the restored 

image. Some additional control over the rate of decay of 

this function away from its centroid is most desirable. By 

_ 2 

making a slight change to Eq. 2-1, replacing [v| by an 
arbitrary weighting function, w(v) , which will be defined 
later, control over the rate of decay of the composite 
imaging system point-spread function becomes possible. 

Thus the fundamental criterion may be expressed as the 
ratio of two quadratic functionals 

" _ 2 _ „ 

/ w(v)g (v) dv 

F i - — — s~rz — i ' (2 - 3) 

/ g (v) dv 
— 00 

The optimum restoration filter would be determined by 
choosing a restoration filter that will minimise the 
numerator of Eq. 2-3 subject to a constraint on the 
denominator, which corresponds to the energy in the 
composite filter. 

In order to avoid two significant shortcomings 
encountered with the classical "inverse" filter, two 
additional constraints must be considered. The first deals 
with the compactness and rate of decay with respect to the 
spatial origin of the point-spread function of the optimum 
restoration filter. Since the point-spread function of the 
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restoration filter based upon the "inverse" filter is 
generally much less compact than the point~ 3 pread function 
of the imaging system, and often is of infinite duration as 
a result of the infinite energy usually associated with 
such a filter, it is advantageous from a data processing 
efficiency and cost viewpoint to reduce the duration of 
this function, thus reducing the amount of data required 
for image processing. In addition since the image to be 
restored usually represents only a portion of the original 
image presented to the imaging system input, significant 
"edge-effect" errors can result. This error is produced by 
truncation of the imaging system output, i.e. having only 
a finite data record available for restoration. It is 
usually greatest at the boundaries or edges of the restored 
image and propagates through the restored image 1 . By 
making the point-spread function of the restoration filter 
more compact about the spatial origin, both problems 
concerning data processing efficiency and "edge-effect" 
errors should be diminished. 

Using an analogous development to that of Eq. 2-3, 

“ _ 2 _ _ 

/ s (v) h r (v) dv 


may be defined, where s (v) 
in purpose to that of w(v) 


oo 2 ' f 

j h (v) dv 

— OO 

is a weighting function similar 
in Eq. 2-3, and h r (v) is the 


i 


See Appendix A 
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restoration filter point-spread function. By using both 
the numerator and denominator of Eg. 2-4 as additional 
constraints, arbitrary control with regard to compactness 
and duration of the optimum restoration filter is possible. 

The second additional constraint necessitated by the 
performance of the "inverse” filter is a constraint on the 
noise power present in the restored image. Because of the 
shape of the spectrum of the "inverse" filter, optimum 
performance is possible only in a noiseless environment. 

The presence of any noise, however small, would in general 
be infinitely amplified by such a filter. The additional 
constraint provides arbitrary control of the signal-to- 
noise ratio of the restored image. This constraint can be 
stated as an expectation over n* (v) , 

F 3 *E{n^(v)} (2-5) 

where n T (v) , the total noise in the restored image, is 
defined in Figure 2-1. 

A final constraint has been found necessary under 
certain conditions to control the shape of the modulation 
transfer function of the composite imaging system to 
guarantee that this function contains no significant 
spurious high frequency peaks which could contribute to 
oscillation near the spatial origin of the composite 
imaging system point-spread function. This occurs 
primarily when the first constraint is removed or is near 
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its unconstrained value, l’his constraint may be defined in 
the spatial frequency domain as, 

® Tf 2 2 

f. = / (“■£■“) |G(D! 

-°° H 

where G(f) is the composite imaging system modulation 
transfer function and f h is a frequency weighting 
parameter. Since the previous constraints were defined in 
the spatial domain, this constraint, as will be shown in 
the next section, can be written as, 

F, * / g"*(v)g(v) dv. (2-6) 

4ir A f 2 ~ m 
H 

2.3 Development of Optimum Filter Function 

Combining the primary criterion defined in Eq. 2-3 
with the constraints of Eq, 2-4 th 2 *ough 2-6, the optimum 
restoration filter point-spread function h r (v) is that 
function which minimizes the functional 

F = / w(v)g (v) dv (2-7) 

—00 

where g(v) is defined in Eq. 2-2, 
subject to the following constraints. 


= / s(v)h 2 (v) dv 

(2-8) 

K 2 ■ E{n 2 (v) } 

(2-9) 


K- / g--(v)g(v) dv. 

J 47t2f n 


(2-10) 


1 
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It should be noted that the constraints related to the energy 
of the composite imaging system, i.e. the denominator of Eq. 
2-3, and the energy of the image restoration system, i.e. the 
denominator of Eq. 2-4, have been dropped since these con- 
straints can always be satisfied by appropriate amplitude 
scaling of the optimum restoration filter defined by the 
constraints of Eq. 2-8 through 2-10. 

Lagrange multipliers and the methods of functional 
analysis will be used to solve Eq. 2-7 subject to the con- 
straints of Eq. 2-8 through 2-10 by forming an augmented 
functional which is the sum of the functional to be mini- 
mized, Eq. 2-7, and the constraint functionals each multi- 
plied by a different Lagrange multiplier [7,8]. 

Thus the augmented functional, I, which must be mini- 
mized with regard to h r (v) becomes. 


I » / w (v)g 2 (v) dv + / s(v)h 2 (v) dv + A 2 E{n£(v)} 


4tt 2 f 2 -» 

ri 


/ g""(v)g(v) dv 


= / w(v) / h (z)h. (v-z) dz / h (u)h. (v-u) dudv 


CO 

+ \- f f s(u)s(u-v)h (v)h (u) dudv 
— 00 


CO 00 

+ A 2 E{/ n{v-z)h r (z) dz / n(v-u)h r (u) du> 

—00 —co ^ 
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oo eo cp 

/ / <S"(t-v) / h r (z)h b (T-z) dz / h r (u)h b (v-u) 

mOO — CP 


dTdvdu. 


( 2 - 11 ) 


Eg. 2-11 way be written in quadratic functional form as 

CO 

I ■ / / A^tiiTraDh (z)h r (u) dzdu 
— 00 


00 

+ X / / A 2 (u # z)h (z)h r (u) 

mOO ■** 


dzdu 


+ X 2 / / A 3 (u f z)h r (z)h r (u) 

• OO 
oo 

+ X 3 / / A 4 <u,z)h r <z)h r <u) 

™oo ^ 


dzdu 

dzdu 


(2-12) 


where A^ (u,z) ,A 2 (u,z) ,A 3 (u,z) , and A 4 (u,z) 

are linear operators defined as [8] 

00 


A x (u,z) = / w(v)h b (v~z)h b (v-u) dv (2-13) 

—00 

A 2 (u,z ) = 3 (u) 6 (u - z) (2-14) 

A 3 (u,z) = E{n(v - z)n{v - u) } 

- R^ (z-u) , f or n ( • ) a stationary 
ergodic random process (2-15) 


and 

A 4 <u,z) 


4u 2 f H 2 -« 


/ / 6'"{t - v) h, (t-z) h, (v - u) dxdv 


■, “> 

/ (v-z) h b (v-u) dv. 


«* 2 V 


(2-16) 
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The stationary points of the augmented functional of 
Eg. 2-12 will be determined by equating the directional 
derivative of this functional to zero. It should be noted 
that for a real functional equating the directional 
derivative to zero is equivalent to setting the gradient to 
zero [8] . The gradient of a quadratic functional 

CO 

I, = / / A(u,z)h(z)h(u) dzdu 

wOO 

is defined as 

oo 

VI, = / [A(u,z) + A'(u,z) ]h(z) dz 

^ oo 

where A'(u,z) is the adjoint of the linear operator A(u,z) . 
The adjoint linear operators of Eq. 2-13 through 2-16 

are 

A£(u,z) = A 1 (u / i‘) (2-17) 

A£(u,z) = A 2 (u,z) (2-18) 

a£(u , z ) * A 3 (u r z) (2-19) 

and a£ (u,z) = a 4 (u,z) . (2-20) 

Setting the gradient of Eq. 2-12 with respect to h r («) 
equal to zero, a homogeneous Fredholm integral equation of 
the second kind, sometimes referred to as a Fredholm 
integral equation of the third kind, is obtained as 

oo 

J [A^VrZ) + X 2 A 3 (v,z) + X 3 A 4 (v,z))h r (z) dz" 
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+ X 1 fi(v)h r {v) = 0. (2-21) 2 

The solution of this equation coupled with the constraint 
equations, Eq. 2-8 through 2-10, would give the required 
point-spread function of the restoration filter, h r (v). 
However because of the numerical difficulties which may 
arise in the general solution of this type of equation and 
in order to more conveniently utilize numerical techniques 
which are readily available, the solution for h r (v) will be 
formulated in the two-dimensional spatial frequency 
domain . 

Eq, 2-12 may be rewritten using inner product notation 
[8], where the inner product of x(t) and y(t) is defined 
as 

00 

(x,y) « / x(t)y* (t) dt 

— oo 

where the asterisk denotes the complex conjugate. 

Thus 

1 = {A l h r' h r 3 + A i (A 2 h r' h r } + x 2 (A 3 h r'V 

* X 3 (A 4 h r' h r K (2-22) 

Eq. 2-22 may also be written as a quadratic functional in 


2 


For derivation see Appendix B. 
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the spatial frequency domain as 

I = + X l (B 2 H r ' H r > + X 2 (B 3 H r ,H r ) 

+ X 3 (B 4 H r/ H r ), (2-23) 

where 

H r (.) « 5'{h r (-)> 

? ( • ) denotes the Fourier transform and , B 3 , and 

are the spatial frequency linear equators which are the 
Fourier transforms of the spatial domain linear operators , 
A^, A 2 , A 3 , and A 4 , respectively. Thus 

B 1 (f,v) = A A (u / J)5i 2,r?5 a; a? (2-24) 

mOQ 

which after substituting Eq. 2-13 into Eq. 2-24 may be 
simplified to 

B-j^v) * H*(?) H b (v)W(f-v) (2-25) 

where H b (») “SKh^*)) 

W(*) = ff{w(-)} 

and H^M*) t ^ ie com Pl ex conjugate of H^(*). 

The adjoint of B^(1T r v) f defined as b£(?,v), may be written 
as 

b£(F,v) = B*(v,f) 

which from Eq. 2-25 becomes 

B£(?,v) « H b (v)H* (1T)W* (v-F) . (2-26) 


Similarly 


B 2 (F,v) = sflf-v) 


(2-27) 
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where 


where 


and 


s(.) - 


B' (f,V) = S*(V-f ) 

(2-28) 

B 3 (f,v) = a> nn (v) (f-v) 

(2-29) 

= ?{R nn ( ’ )} 


B'(S,v) = $ nn (f)6(v-f) 

(2-30) 

B 4 (f,v) - ~ |i^(v)| 2 6(f-v) 
f H 

(2-31) 

Bj(f,v) = I (f ) j 2 S<V-f). 

f H 

(2-32) 


The gradient of the quadratic functional of Eq. 2-23 
becomes 

VI = (B 1+ Bj)H r + X l( B 2 +B')H r + X 2 (B 3+ B')H r 

H X 3 (B 4 +B')H r (2-33) 

which upon expanding the linear operator notation of Eq. 

2-33 becomes 

CO 

VI - / [B x (f,v) + Bj{f,v) + X 1 {B 2 (f,V) + B'(f,v)} 

wCO 

+ X 2 (B 3 (f r v) + B^{f r V)> 

+ X 3 {B 4 (f r v) + Bj(f,v)}]H r (v) dv = 0. 

(2-34) 


3 


See Appendix C for a derivation of these spatial frequency 
linear operators. 
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Substituting Eg. 2-25 through 2-32 into Eq. 2-34, 

vi «= / [h. ffm (vHw(f-v) + w* (v-f ) } 

— Cq 

+ A 1 {S(f-v) + S*{v-£)}]H r (v)dv + 2[A 2 § nn (f) 

+ X 3 (|-) 2 1^(1) i^H^f) = 0. (2-35) 

Eq. 2-35 represents the general expression for the 
gradient of Eg. 2-23 with respect to H (•), which then 
combined with the constraint equations, completely specified 
the spatial frequency spectrum of the restoration filter. 

The constraint equations, Eq, 2-8 through 2-10, may be 
rewritten in the spatial frequency domain as 

K x = (B 2 H r ,H r ) (2-36) 

which after substituting Eq. 2-27 into 2-36 becomes 

K 1 = / /S(f-v)H r (v)H r *(f)dvdf, (2-37) 

“CO ” 

K 2 = (B 3 H r ,H r ) (2-38) 

after substituting Eq. 2-29 into 2-38 becomes 
00 

K 2 = / § nn (?) |H r (f) | z df (2-39) 

— oo 

and 

K 2 = (B 4 H r ,H r ) (2-40) 

after substituting Eq. 2-31 into 2-40 becomes 
K, = / C|— ) 2 | H. (f > | 2 ! H (f)| 2 df. 

J -to H 


(2-41) 
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Before Eq. 2-35 can be reduced to a form more suitable 
for the evaluation of the spatial frequency spectrum of the 
restoration filter, the penalty weighting functions w{v) 
and s (v) , introduced in Eq. 2-3 and 2-4 respectively, must 
be further examined. Since w(v) is designed to influence 
the solution of h (v) so that the composite imaging system 
point-spread function, g(v>, is duration limited, a 
possible choice for w(v) would be 


w (v) =1 f or v^ <_ v 


— v 2 


= otherwise , 

However such a choice for w(v) would lead to analytical 
difficulties in Eq. 2-35, since the Fourier transform of 
such a function does not exist. Thus the expression for 
w(v) must be chosen in such a manner that it allows enough 
flexibility to arbitrarily control the duration or radius 
of gyration as well as the rate of decay of g(v) and, in 
addition, to have a Fourier transform [25] . 

One function for w(v) which satisfies the previous 
requirements, written in terras of one variable, is [15] 

2v - v 1 - v 9 2k 

w(v) - [ — j ■■■■-] + e (2-42) 

v 2 V 1 


for 0 < c < 1 

and k a positive integer. This function is shown in 
Figure 2-2 . For convenience s (v) will also be described by 
the same type of function. 
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The following analysis is based upon a rectangular 
coordinate system. Eq. 2-35 may be rewritten in terms of 
the x- and y- components of f and v. 


// [H. * (f ,f )H. (v ,v ){W(f -v ,f -v ) + W*(v -f ,v -f )} 
L ~b ' x' y d x y x x y y x x y y 


+ A,{S(f -v ,f -v ) + S (v -f ,v -f ) } ] H (v r v )dv dv 
1 x x y y x x y y r x y x y 


2 . xs 2 


f A+f. 


+ 2 [X„t (f ,£ ) + A - {— >X_)|E(f ,f )I 2 ]H (f ,f ) « 0 

2 nn' x y 3 -2 1 n x' y 1 J r x y 


f 2 

H 


y 

(2-43) 


For convenience in handling the analysis with respect to a 
rectangular coordinate system both w(v) and s(v) will be 
defined as the product of their x- and y- components, from 
Eq. 2-42, 

w(v) = w (x) w (y) (2-44a) 

x y 


2x-x . — x 0 2k 2y-y _-y _ 2k 

= r{ y- 1 w2 > + c ] [{ wy + c ] 

x w2" x wl x Y w2” Y wl Y 

(2-44b) 

where 0 < c <1 

X 

0 < c <1 

y 

and for k and k w positive integers, 
wx wy 

Similarly 


s(v) = s (x) s (y) 
x y 


(2-45ar) 


= [{- 


2x-x ,-x „ 2k 
si s2j sx 


x s2" x sl 


■2y-y , -y „ 2k 


d 3 gZ sl_ s2 } S y + d 3 

x v_o-y_-i v 


s2 J sl 


y 

(2-45b) 


I 
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where 0 < d x < 1 

0 < a y < 1 

and for k and k positive integers, 
sx sy 

Choosing k wt =» 1 = , the Fourier transform of 

Eg. 2- 4 4b becomes 


W(f ,f ) = W (f )W (f ) 

v x' y x x y' y 


where 


W (f ) = 
x x 




( 2 - 46 ) 

, j 2 (x .+x J 

[- Ji_ ($■'•'(£ ) — *2- 6 •* (f ) 

2 x' 7T x' 

7 T 


and 


+ t(x wl +x w2> 2 + (x w2- x wl )2c x ) S{£ x ,] (2-«) 


1 i 52 (y ,+y 2 ) 

W (f ) = t- — 6 "" (f ) S'(f ) 

y y <^2-y w i> 2 » a y y 


{(y wl +y w 2> 2 + ( 2 - 48 ) 


w 2 J wl' y J 


A similar expression results from taking the Fourier 
transform of Eg. 2-45b, 

S < f x'V “ S *‘WV ( 2 - 49 ) 


where 


1 i ] 2 (x _+x 

S (f ) = ± [- — fi'-'ff ) — 6 (f ) 

x x (X -x ) 2 IT 2 X ^ x 

lX s 2 X sl J 17 


+ {<X S 1 +X S 2 )2 + <x S 2 ‘ x sl ) 2 d x ){i(, x 11 ( 2 - 50 ) 
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and 


i i j2(y_ 1 +y_o) 

s (f ) = I - ~ 6 ■"’ ( f ) s ;= <$'(f ) 

Y y <y S 2-y s i' 2 * 2 y y 


+ { ( y s i + y S 2> 2 + ( y S 2-y B i ,2a y }4( V l (2 ' 51) 


If it is assumed, as would usually be the case, that the 
penalty functions, w(x,y) and s(x,y) are centered about 
the origin, then 

x wl “ ■ x w2 " x w 


y w i " - y W 2 ■ * 


w 


x sl = - X s2 - x s 


*sl ~ ^s2 = 


(2-52) 


Substituting Eq. 2-52 into 2-47, 48, 50, and 51, 


W (f ) = - 

x x 


4tt x 


-±_ S"tf x ) + c x S(f x > (2-53) 


w 


Wff) = - (f ) + c 6 (f ) (2-54) 

y y 4-?r 2 y 2 y y y 


w 


S (f ) = - — 

x x 47f2x s 2 6 "" (f X ) + d x 6 (f x } 


(2-55) 


and 


S (f ) = 

y y 


8 " (f ) + d 6(d) 

4ir 2 y 2 ^ y y 


2 


(2-56) 
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Substituting Eg. 2-53 through 2-56 into Eg. 2-43 


gb <W 

8tt 4 x 2 y 2 
w J w 


3f 2 3f 2 
x y 


vvyvvy 


- 4ir 2 y 2 c 
J w y 

3 2 

H. (f ,f ) H (f ,f ) 
a x' y r x y 

3f 2 

X 

- 4tt 2 x 2 c 
w X 

3 2 

H. {£ ,f )H (f ,f ) 
o x' y' r x y 

3f 2 

y 

+ 8tt 4 x 2 y 2 
w J w 

C-CHL (f ,f ) H (f ,f )} 

x y b x y r x v 


l l f 3" 


8tt 4 x 2 y 2 3f 2 0f 2 ’ V f x' f y } 
s J s x y 


- 4ir 2 y 2 d 


H (f ,£ ) 


s y 3f 2 r x' y 


x 

s 2 


- 4it 2 x 2 d — — H (f ,f ) 

s x 3f ^2 r x' y J 


+ 8Tr‘ t x„ 2 y 2 d d H (f ,f ) } + 2[X„4„ {f ,f ) 
s - r s x y r x y 1 2 nn x y J 


(f 2 +f 2 ) 

+ A- -* X— 

3 f 2 

h 


IVVV^VW ■ °' 


which when expanded in terms of H r (f x ,fy), becomes 


[A(f x'VV f x' f v> +Btf x' f v )] 


3 4 H (f. ,£ ) 
r x y 

3f 2 3f 2 
x y 


3K (f ,f } 3 3 B (f ,f ) 

_bj x' y' r x' y 


3f 


3f 2 3f 

^ y 


(2-57) 


+ 2A(f ,f ) 

x y 
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+ 2A(f x , 

+ 4A(f x , 

+ [A{f x , 

+ G ( f x t f 

+ D(f ,f 
x 

+ [2A(f 

x 

+ 2E(f x ,: 
+ [2A(f 

A 

+ 2D(f x ,: 

+ [C( V 

+ E(f x ,f 


f ) 

y 


V 


f ) 

y 


aH b ( W 


9f 
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(2-58) 


where 


A ^ f x' f y^ 

ts 

=b < f x'V 

(2—59a) 


- 

*A* { VV 

8^x s *y s 2 

C 2 — 5 9 b ) 

C{f ,f ) 
x y 

- 

(f z +f z ) 

IO xPy + 2X 3 * f / 1 1 W 

E 

f ) r 
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+ Wy + 2 ‘ 2 y £ S 'V 

(2-59c , 

D (f ,f ) 
x' y' 


<W°V 

(2-59d) 


2tt z x 2 
w 

E(f , f ) 
x y 


-HU * ( f f f ) c 
& x r y' X 

( 2-59e ) 


2ir z y 2 
J w 

F(f ,f ) 
x y 


-Aid 

2tt z x 2 
s 

(2~59f ) 

GtE x'V 


X 

nd 

. i— 1 
r< 

1 

(2~59g ) 


2ir 2 y 2 
s 

H(f x ,£ y ) 

- 

c c III. (f ,f ) I 2 + A, a d . 
x*y 1 b x y 1 1 x y 

(2-59h } 


Similarly in terms of the rectangular spatial frequency 
coordinate system the constraint equations, Eq. 2-37, 39 , 


and 41 become, 
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K„ = 


K„ = 


//// S(f -v r f -v ) H {v ,v )H (f ,f )df df dv dv 
- oo ; x x Y y r x' y r x y x y x y 


// $ (f ,f ) j H {f ,f ) | 2 df df 
nn x y 1 r x r y ' x y 


(2-60) 

(2-61) 


K- = 


°o f 2 + f 2 
// ( x 


2-) Is (f ,f ) I 2 I H {f ,f ) I 2 df df . 

d x r y 1 1 r x y 1 x y 

r H (2-62) 


Thus Eq. 2-53 and 59 in conjunction with the constraint 
equations, Eq. 2-60 through 2-62, specify the general form 
of the optimum restoration filter in terms of the rectangular 
spatial frequency coordinate system. 

2.4 Separable Optimum Restoration Filter 

The general system of equations defined by Eq. 2-58 
through 2-62 represents a formidable solution problem, 
which in general would be possible only by a two-dimensional 
numerical solution technique requiring a computer system 
having very large data storage capabilities. However for a 
certain class of aperture functions, which represent a large 
class of physically realizable apertures, an enormous 
reduction in the computational complexity of the problem is 
possible. This class comprises those apertures which may be 
modelled as separable apertures, where it is assumed that 


^x ' f y ^ ^bx^x^y (f y } 


§ (f ,f ) = $ (f )4> (f ) 
nn x ' y nnx x nny y 


H (f ,f ) - H (f )H (f ). 
r x' y' rxxryy 


(2-63) 

(2-64) 


(2-65) 
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With the assumptions of Eq. 2-63 through 65, the solution 

of Eq.2-58 through 62 can be considerably simplified by 

use of the method of separation of variables. Instead of 

substituting Eq. 2-63 through 65 into Eq. 2-58 and 

separating Eq. 2-58 into two differential equations, one 

a function of f and the other a function of f , a somewhat 
x y 

more fundamental approach will be used. 

Taking the inverse two-dimensional Fourier transform 
of Eq. 2-63, 

h b (x,y) = h bx (x)h by (y) (2-66) 

Similarly Eq. 2-64 and 2-65 become respectively 

R im (x ' y> = E imx (x)R nny (y) < 2 -6’> 

and 

h r (x,y) = h rx (x)h ry (y) • (2-68) 

Writing Eq. 2-7 in terms of the two spatial dimensions, 

CO 

F = ff w(x,y) g 2 (x,y) dxdy (2-69) 

— CO 

where from Eq. 2-2, 

g(x,y) « h r (x,y)**h b (x,y). (2-70) 

Substituting Eq. 2-66 and 2-68 into 2-70, 


g(x,y) = [h rj£ (x)h ry (y) ]**[h bx (x)h by (y) ] 
where 

g K (x) “ h rx <x) * h bx (x > 


g x (x)g y (y) 
(2-71) 


(2-72) 
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and 



g y (Y) = h ry (y) * h by (y). 

(2-73) 

Substituting 

Eg. 2-44a and 2-71 into 2-69, 



It 

(2-74a) 

where 

00 



F x = / w x (x)g^(x;dx 
—00 

(2-74b ) 


CO 

F y = / w y (y) g y (y) dy . 

(2-74e) 

Substituting Eg. 2-45a and 2-68 into 2-8, 



K. = K. K n 
1 lx ly 

(2-75a) 

where 

CO 

K lx - / s x U)h 2 rx ( X )dx 

“00 

(2-75b ) 

and 

K iy - j s y (Y,hZ ry (Y)dY - 

(2-75c ) 

Substituting 

Eq. 2-64 and 2-65 into Eg. 2-39, 



K 2 - K 2x K 2y 

(2-76a) 

where 

CO 

K- = / $ (f ) |H (f ) 1 2 df 

2x J nnx x 1 rx x 1 x 

— CO 

(2-76b) 


and 

co 

K 0 = / § (f ) |H (f ) | 2 df , 

2y nny y 1 ry y 1 y 


(2-76c) 
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and from Eq. 2-15 


K 2x = E te 2 Tx (^)> 


(2-76d) 


and 


K 2y = E{n2 Ty < y >> 


(2-76e) 


And finally writing Eq. 2-6 in terms of the two spatial 
dimensions, 

K — - // g" (x,y)g(x,y)dxdy (2-77) 

3 4TT 2 f 2 -co 

£1 

and after substituting: Eq. 2-71 into 2-77 , 


K f g ■"* (x) g (x)dx /g "^(y)g (y)<3y. 

3 4ir 2 f e 2 -» x x -» y y 


(2-78) 


Using the results of Eq. 2-12, 73, 16, 20, 31, 32, and 41, 


where 


K 3 - K 3xSy 


« f 


(2-79a ) 


K 3x - ^x< f xH^H rx C f x )| 2 df x (2-79b) 

-co HX 


and 


00 f 


K 


3y 


/ CjjM l^ y (* y )l I l B ry ( Vl' M y (2 - 79c) 


Hy 


where 


f 2 = — 4ir 2 f 2 f 2 

H Hx Hy 


(2-79d) 
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Thus the problem of determining the optimum 
restoration filter point-spread function, h (x,y) reduces 
to finding the h (x) that will minimize Eq. 2-74b subject 

JL Jim 

to the constraints of Eq. 2-75b, 76b, and 79b, and to 
finding the h ry (y) that will minimize Eq. 2-74c subject to 
the constraints of Eq. 2-75c, 76e, and 79c. The original 
two-dimensional restoration problem reduces to two one- 
dimensional processes which have similar equations. From 
the preceeding vector notational analysis used in Eq. 2-7 
through 2-41, 2-44, and 2-45, the system of equations 
necessary to solve for h^(x) and h^Cy) may be formulated. 

To solve the h (x) , the augmented quadratic functional 
of the form of Eq. 2-11 determined by Eq. 2-74b, 75b, 76b, 
and 79b becomes, 

DO 00 

T x = / w x (x)g 2 (x)dx + X / s x (x)h rx 2 (x)dx 

* —CO ‘ —00 


+ *2x E{n z Tx (x)} - - .* 3x . / (xjg^txldx 


«* 2 V - X 


(2-80) 

which may be written as a quadratic functional in the 
spatial frequency domain from Eq. 2-23 as 


T x - < E lx H rx' H rx> + X lx (B 2x H rx' H rx> 

+ i 2x< B 3x H rx' H rx> + X 3x (B 4x H rx' H rx> 


(2-B1) 


where 


B lx t£ x'V - '‘bx t£ x> H bx ( W f x-V 


(2-32) 
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from Eq. 
from Eq. 
from Eq« 

from Eq. 

from Eq. 

from Eq. 
and from 


The 
Eq. 2-33 


2-26, 


H bx < f x> W x ( V f x J 

2-27, 


(2-83) 


B 2x (f x' V x> 

2-28, 


S (f -V ) 
x x x 


(2-84) 


B 2x (f x'V 


S (v -f ) 
x v x x 


(2-85) 


2-29, 


B (f #V ) = (v ) 6 (f — v ) 

a 3x' x x^ nnx x' x x' 


( 2 - 86 ) 


2-30, 


B 3x' (f x’V “ W f x> “W 


(2-87) 


2-31. 


B 4x (f x'V - <f3 S( V v x> (2- 8S > 


Hx 


Eq. 2-32, 

B 4x' (f x'V " IV (f x)l 2 «<V f x> 


(2-89) 


gradient of Eq. 2-81 may be written in the form of 


7I x = (E lx +B lx' )H rx + A lx (B 2x +B 2x' )H rx 


+ *2x (B 2x +B 3x ,H rx + A 3x (B 4x +B 4x )H rx 


(2-90) 
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which may be expanded to the form of Eq. 2-35, 

00 ^ 

VI = f [H. (f ) H, (x){W (f -V ) + W (v -f }} 
x L ox' x' £x' 7 X v X X X X x' 


+ WWV + S X ( ' , x- f x !}IH r X t ' ,)av X 

r Hx (2-91) 

Substituting Eq. 2-53 and 55 into Eq. 2-91, the following 
differential equation results 

H ™ (f x > + i 2 <x< f x> B £*« £ x’* B * /I l H bx‘Vl 4 x s ,+X lx 1 S,’ 1 }B ^ t£ x> 

+ tx s 2 'C l£ x >^ (£ x > - 4 ” 2x w i{c x w 3 X ( y-> 2 'l^ x (f x ) l 2 

Hx 

+ X 2xWV +d x }] / ll^«x>IV +X lA a3, ^*< £ x> 


= 0- 


( 2— 9 2 ) 4 


From Eq. 2-36, 2-55, and 2-75b, 


K 


i-/ IK— (f )H _(£ ) + H'l (f )H . (f )]df 


lx . 2 2 J rrx ' x' rrx' x' nx' x' nx' x‘ 

4tt x_ - m 


+ d / | H (f ) I 2 df 
x J 1 rx x 1 x 
—00 


(2-93) 


4 For derivation, see Appendix D. 


where 


H rx^ f x^ H rrx^ f x^ + ^ H rix^ f x^ (2-94) s 

and where H (f ) is the real part and H . {£ ) is the 
rrx x ^ rxx x 

imaginary part of H (f ). 

xTX X 

Restating Eg. 2-76b, 

K 2* ' /• nn X (f X >l H r X « X ’l 2af X < 2 -™b) 

— DO 

and Eg. 2-79b, 

K 3x " l^ X < f X >l a |H rx tf x )| 2 df x . (2 _ 79b) 

-W HX 


Thus the simultaneous solution of the differential equation 
Eg. 2-92 f and the constraint equations, Eq. 2-93, 2-~76b, 
and 2-79b, specify the form of the x- component of the 
spatial frequency transform, or equivalently the point- 
spread function, of the restoration filter. 

In a similar manner it is possible to solve for h^ {y ) 
by forming an augmented quadratic functional of the form of 
Eq. 2-11 determined by Eq. 2-74c, 75c, 76e, and 79c, 

x y *• _/ w y (y)g y 2 ty)dy + * ly _J s y (y)h ry 2 (y )dy 


+ X 2y E ^ nZ Ty <Y)> " 


Jl. 


4ir 2 f 2 

x Hy 


/ g v '"'(y)g v (y)ay. 

» y y (2-95) 
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For derivation, see Appendix E. 



By following an analogous procedure to that used for 


determining h (x) in Eg. 2-81 through 2-94, the equations 

IT A 


which specify h r y(y) may be formulated. Only the results 


will be stated since the derivation of the equations for 


h (y) is identical in procedure to that given for h (x) 
ry rx 


with the appropriate change in variables form x- to y~ 
dependency . 

The differential equation specifying the form of 


H ry (fy) becomes. 


n^<f y ) + {z^vvv^/i'vv iv + H y v ], vy 


+ ty s at ^y (f y ,E £y (f y > ' "VV 1 


+ $ (f )+d }] /I I h (f ) 1 2 y 2 +X-, y 2 ] } h (f ) 

2y nny y y j/1 - , r>y'y' |jr s . ly^w ry y 


= 0, 

while the constraint equations become, 

CO 

1 


(2-96) 


K 


iy 


f 2 v 2 
S 


4ir z y ‘ -oo 


/ [E" (f ) H (f ) + H"r (f ) H . (f )]df 

J TrTTXT X7 ' T'T'W XT V XT * ri \7 XT 1 


+ d / t h (f ) 2 df 

y r y y y 


rry y rry y riy y ny y y 

(2-97) 


K 0 — / $ (f ) I H (f)[ 2 df 
2y _£ nny y 1 ry y 1 y 


(2-76c) 


and 


oo f 


K 3y=J <S^> IVV^VV^ 


(2-79c) 



Thus the simultaneous solution of Eq. 2-96, 2-97, 2-76c, 
and 2-79c will specify the y- component of the restoration 
filter. 


2.5 Radially Symmetric Restoration Filter 


Probably the most common type of aperture, because of 
physical ease in construction, is the radially symmetric 
aperture. For this case it is assumed that 


VW 




(2-98) 


0 (f ,f ) = 

nn x y 


$ (f ) 
nnr r' 


(2-99) 


H (f ,f ) 
r x y 


E rr (f r> 


( 2 - 100 ) 


where 


f 2 +f 2 

x y 


( 2 - 101 ) 


The solution to Eq. 2-58 is analogous to the solution 
for h rx (x) where all variables dependent upon x are 
replaced by corresponding r dependent variables. From 
Eq. 2-44a, 


w(v) = w r (r) 


(2-102a) 


2r-r ,-r , 2k 
_ . wl w2 -j wr 

r w2 wl 


(2-102b) 


for k wr a positive integer and 

r 2 = x 2 + y 2 * 


(2-102c) 
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Prom Eq. 2-45a, 

s(v) = s r (r ) 


2r~r _-r „ 2k 
_ j. si s2 -] sr 

r s2*~ r sl 


(2-103a) 

(2-103b) 


for k gr a positive integer. 

Choosing 

k * 1 * k and 
wr sr 

r wl = " r w2 * r w (2-104) 

si s2 ^s' 

the Fourier transforms of Eq. 2-10 2b and 103b become 


W (f ) 
r r ' 



s~(f r ) 


(2-105) 


S (f ) = J™ 6 '"'(f ). (2-106) 

r r 4iT 2 r 2 r 

s 

To solve for h (r ) , the augmented quadratic 
functional of the form of Eq. 2-11 which must be minimized 
with respect to k rr (r) becomes, 

CO co 

J r ” / w r (r)g r 2 (r)dr + X^ / s r (r )h rr 2 (r )dr 




*2r ={n 2 Tr (r)} - 




f (r)g r (r)dr . 

(2-107) 


By following an analogous procedure to that used for 
determining h„„(x) in Eq. 2-81 through 2-94, the equations 


I 
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which specify h^fr) may be formulated. Again only the 
results will be stated since the derivation of the 
equations for h r ^(r) is identical in procedure to that 
given for h (x) with the appropriate charge in variables. 

The differential equation specifying the form of 
H rr ( f r ) becomes. 


H « (f r> + {2 <r (f r>^r tf r> r s ^ [ l B br (£ r ) I 2r S 2+X ir r w 2,lH i r (f r> 


f r 2 


+ fr s 2 ' H br (f r , ^< f r> " 4 *V {X 3r ( 'VVl 

nr 


+ A 2r*nnr tf r ))] ' 1 1 V (f r ) I V+hrV 1 > H rr (f r> 


= 0 


(2-108) 


while the constraint equations become. 


K 


lr 


CO 

/ f H ^r (f r )H rr< f r> + H ^r < f r > H rir <«r> 


4iT 2 r “ 
s 


K 0 = / 0 (f } lH (f ) | 2 df 

9 t ^ nnr ' r* 1 ' irir r* 1 t 


and 


K 3r = i ( J^ ,2 | B br (f r>l 2 l H rr tf r>l 2df 1 


nr r nr r 

(2-109 ) 

( 2 - 110 ) 
( 2 - 111 ) 


Thus the simultaneous solution of Eq. 2-108 through 2-111 
will give the spatial frequency spectrum of the optimum 
restoration filter. 
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2.6 Additional Shape Control of Composite Imaging Function 
by a Sub-optimal Iterative Technique 

The third constraint, originally defined by Eq. 2-10, 

was chosen to provide control over the composite imaging 

function in such a way that the magnitude of secondary 

lobes in the immediate vicinity of the primary lobe could 

be reduced. These secondary lobes lead to undesirable 

ghosts in the restored image. Such a constraint is 

necessary because the weighting function, w(v) defined by 

Eq. 2-42, of the fundamental design criterion as expressed 

by Eq. 2-7 has almost no significant effect on g(v) in the 

immediate vicinity of the spatial origin. Thus, without 

the third constraint, significant secondary lobes in this 

region of the spatial domain could exist. 

In order to provide additional control over the shape 

of the composite imaging function in the vicinity of the 

spatial origin, which is not feasible by the third 

constraint alone, an iterative technique is utilized. To 

implement the technique, the fundamental design criterion, 

defined by Eq. 2-7, is modified to be, 

00 

F = / W (v) [g (v)-m(v ) ] 2 dv , (2-112) 

— CQ 

where the criterion function m(v) will be defined shortly. 

The procedure for this technique is as follows. 
Initially, m(v) is set equal to zero; Eq. 2-112 thus 



I 


J 
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reduces to the original fundamental criterion as defined 
by Eq. 2-7. A solution for the restoration filter based 
upon minimizing Eq. 2-112 subject to the constraints of 
Eq. 2-8 through 2-10 is obtained. The existence of a 
unique solution to this system of equations will be 
discussed in the next chapter. It should be noted that the 
resulting composite imaging point-spread function, g(v), 
is optimal with regard to having a minimum radius of 
gyration. However g{v) may still have some undesirable 
side lobes in the immediate vicinity of the spatial origin 
which were not sufficiently suppressed by the third 
constraint. In such a case, m(v) is defined to be equal 
to the optimum g(v) in the region about the origin 
containing the primary response lobe but excluding any 
secondary side lobes, and equal to an appropriately chosen 
function of the optimum g(v) outside of this region so as 
to further constrain against the side lobes. 

After choosing an appropriate m(v), the first iteration 
is made by minimizing Eq.. 2-112, subject to the constraints 
of Eq. 2-8 through 2-10. The resulting composite imaging 
function, although no longer optimum with regard to having 
a minimum radius of gyration, will not deviate significantly 
from the optimum in this respect because of the way in 
which m(v) is chosen and-, more importantly, will have 
significantly lower secondary lobe responses due to the 
weighting effect introduced by m(v). If necessary, a second 
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iteration may be made by choosing an appropriate m{v) based 
upon the g(v) of the first iteration, and so on. There 
appears to be an almost unlimited number of functions which 
may be chosen for defining the criterion function including 
the absolute value, and exponential, Gaussian, and poly- 
nomial approximations; however the effect of any of these 
upon the iterated composite imaging point-spread function 
should be similar, because of the way in which m(v) is 
chosen . 

An example of one possible one-dimensional criterion 
function is shown in Figure 2-3. For this case, m(v) is 
equal to the optimum g(v) in the region of the primary lobe 
about the origin where g(v) is positive. Outside of this 
region, m(v) is chosen to be the absolute value of the 
optimum g(v). A comparison of the g(v) functions, both 
before and after iteration, would reveal that the secondary 
side lobes have been heavily constrained, while the general 
shape of the primary lobe has remained essentially the same. 
Thus this technique provides a very powerful means for 
making minor adjustments to the shape of the composite 
system point-spread function, g{v). Although g(v) could be 
altered in the same manner as that provided with this 
iterative technique by the introduction of additional 
constraint equations, the complexity of the functions 
necessary to produce the same control over g(v) makes this 
alternative less attractive. 
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The following analysis is based upon the assumption 
of either a separable or a radically symmetric aperture 
since these are the most common apertures and are the forms 
considered in the remainder of this work. The general 
analysis, based upon a non-separable two-dimensional 
aperture is straightforward and analogous to that used in 
the development of Eg. 2-58 through 2-62 except that 
Eq. 2-58 is non-homogenous . 

For a rectangularly separable aperture, Eq. 2-112 will 
be redefined in an analogous manner to Eq. 2-74 as 

F = F x F y (2-113a) 

where 

CO 

F v = / w (x) Eg (x)-m (x)] 2 dx+p (2-113b) 

— CO 

and 

CO 

Fy = / w y {y) [g y (y)-m y (y)] 2 dy+p y (2-113c) 

where 

0 < P x ,P y « If (2~113d) 

The constraint equations remain as defined by Eq. 2-75, 76, 
and 79. Thus for a given ni (x), h <x) must be determined 

4 A Ju A 

by minimizing Eq. 2-113b subject to the constraints of 
Eq. 2-75b, 76d, and 79b and for a given m y (y), h ry (y) must 
be determined by minimizing Eq. 2-113c subject to the 
constraints of Eq. 2-75c, 76e, and 79c. 

To solve for h (x), an augmented quadratic functional 

1TX 

formed from Eq. 2- 11 3b and Eq. 2- 7 5b, 76d, and 79b can be 


written 


CO 


l x = . / w x (x) [g x (x)-m x (x)] 2 ax+ p x 
*“bo 


+ A lx / s x U)h rx 2 (x)dx 

— co 


A 00 

+ * 2x E{n 2 .< x )> - -■ , , / g''(x)g U)dx. 

r3> 4lT f Hx (2-114) 

Expanding the first integral of Eg. 2-114 gives 


I x * / w x tx ^ x 2 (x)dx -2 / w x (x)g x (x)m x (x)dx 

— CO mCO 


+ / w x (x)ra x 2 (x)dx + P x + X lx / s x (x)h rx 2 (x)dx 


+ X, E{n 2 (x)} - 
2x Tx 


3x 

2 7“T L y x 


4TT i! f £ — 00 
Hx 


/ g""(x)g (x)dx. 


(2-115) 


Since the functional form of Eq. 2-115 is basically 
similar to that considered in Eq. 2-80 with the exception 
of the second, third, and fourth terms, only these terms 
will now be considered. The second term in Eq. 2-115 may 
be expanded into operator form as 

co 

I In = — 2 / w x (x)g x (x)m x (x)dx 

— CO 

CO CO 

= - 2 / w x (x) / h rx (T)h bx (x-T)dT m x (x)dx 

—00 —00 ” 

= l A 5x (T)h rx (T> dT ' 


(2-116) 


where 


A Sx (t) = 


- 2 / w x (x)m x (x)h bx (x-T)dx. (2-117) 


Alternately Eg. 2-116 may be written in the frequency 
domain as 


I = / B_ (f ) H (f )df 

1X1 * h v ' t* v v 


r- J- J U V J- i 

5x x' rx x x 


(2-118) 


where 


B sx' f x> “ - 2 E bx' £ x> / ™x‘V»x*‘ f x + V d ' , ! 


and 


M x ( • ) = ? {m x { • ) ) • 


(2-119 ) £ 


Taking the gradient of the linear functional of 
Eq. 2-118 with respect to H (f ) gives the conjugate 

X A A 

kernel B* x (f x ). 


VI = B_ (f )• 
m 5x x 


(2-120) 


Substituting Eq. 2-53 into Eq. 2-119, 

B 5x<*x> - - 2I W f x> l f- TZ7TT *"<*,1 + °x 5 <V> 

— 00 47T X 

w 

M (f +v }dv 

XXX X 


- I W £ x )[ rx t 


d 2 * * 

— — M (f ) - 2c M (f )). 

2tt 2 x 2 df 2 x x x x x /J 


w x 


(2-121) 


6 


See Appendix F for derivation 
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The third and fourth terms of Eq. 2-115 may be 
neglected since they are not functions of h (*)a and thus 

i A 

the gradient of these terms with respect to h ( • ) is zero. 
The gradient of Eq. 2-115 may be written in the form of 
Eq. 2-90 as, 


7I x = (B lx +B lx' )H r X + A l x (E 2 X +B 2 X ' )H r X 


+ X„ (B, +B, **)H + X- (B. +B. ")H 

2x 3x 3x rx 3x' 4x 4x rx 


+ B 


5x 


( 2 - 122 ) 


which, when expanded, is similar in form to Eq. 2-92, 
except that the resulting differential equation is 
nonhomogeneous , 


B r X " (f x> + 2 ^ X (f X > B b X ' tf X )x s 2 B r X ' (t X >/'|H bx (£ x )| ! x 

- hxV ! + x s 2 t4< f x »^ x " (f x ) 

- 4 ’ r V {[ ° X +A 3 X ( W )2l l i tx (f X >l 2 


2 ,, 2 
S 


+ x 0 $ (f )+a } ] h (f )/[|e (f )| 2 x 

2x nnx x x rx x ' 1 ox x 1 


+ X, x 2 ] = 2x 2 x 2 x 2 B c (f ) / 
lx w s w 5x x 


2„ 2 
S 


[Jh (f ) I 2 x 2 +X x 2 ] . 
1 bx x 1 s lx w 


(2-123) 
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Thus Eq. 2-123 in conjunction with constraint Eq. 2-93, 
76b, and 79b specify the shape of the x- component of the 
frequency spectrum of the restoration filter. A similar 
equation set may be written for the y- axis frequency 
spectrum component of the restoration by using Eq. 2-97, 

76c, and 79c, and replacing the x- variables in Eq. 2-123 
with appropriate y- variables. 

For a radially symmetric aperture, Eq. 2-112 may be 
defined as 

CO 

F = / w r (r) [g r (r)-m r (r)] 2 dr, (2-124) 


where w r (r) is as previously defined in Eq. 2-102. The 
augmented functional to be minimized is similar in form to 
Eq. 2-107, 

03 CO 

X * / w r (r)g r 2 (r)dr - 2 / w r (r )g r (r)m (r )dr 

—CO —CO 

00 CO 

+ f w r (r)m r 2 (r)dr + X 1r j s r (r )h rr 2 (r)dr 

X. 


lr ' _ r ' ’"rr 

— CO 


+ *2r Et » , Tr (r,} - - 

4ir z f;, -« 


/ 9'' (r)g r (r)dr. 

Hr ” (2-125) 


Following a similar development to that for Eq. 2-115, the 
resulting differential equation in terms of H rr ( f r ) 


becomes. 


I 


I 


I 
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H"(f r ) + t2< r (f r )^ r (* r ^ s Vl|^ r (f r )l 2 rs 2 +» lr r w ! 


+ {r s 2f 4 (f r>^ (f r> - 4 ’ r V t> 3r ( ir- )2 | B br (f r>l 

Hr 

+ l aW f r»l / t|^ r ' f r>l 2 r s 2 W lr V I}B rr«r> 

= 2* 2r s V B ^tf r ) / H Ii hr< f r , l 2r s 2+,l lr r w 2] 

d 2 * 

B 5r<V ■= zfrf- ITT W’ 


where 


2ir 4 r 


w 


df 


The constraint equations remain as defined by Eq. 5 
through 2-111. 


2 

(2-126 ) 
(2-127) 
!-109 
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2.7 Classical “Inverse 11 Filter 

Intuitively one would expect that if the constraints 
are removed and if the parameters of the weighting 
function, w(v), are appropriately chosen so that only the 
origin is unweighted thus producing an impulse for the 
composite system point-spread function, then the resulting 
filter should have the form of the classical "inverse" 
filter. 

For notational convenience, the following analysis 
will be made in terms of the radially symmetric aperture. 
From Eq. 2-102 for k wr =l and for a symmetric weighting 
function as defined by Eq. 2-104, 

w (r) - <f-) 2 . (2-128) 

r x w 

It is clear that as r approaches zero, only the origin 

W 

remains unweighted in Eq. 2-107. Prom Eq. 2-108 which 
specifies the optimum restoration filter in the frequency 
domain, if approaches zero and all the constraints are 
relaxed, the resulting differential equation becomes 


2E? (f ) 

h j;<v * 


+ 


H^(f ) H (f ) 
br r rr r 




= 0 


or 


H br (f r ,H ^ (f r ) + 2I £r t£ r )H ^<V + " 0 

(2-129) 

If it is assumed that 


H (f ) = k H, (f ) , 
rr r hr r ' 


(2-130) 
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where k is an arbitrary constant, then it can be shown that 
Eq. 2-129 is satisfied. Since there exists a unique solution 
for a given differential equation of the form of Eq. 2-129 f 
then Eq. 2-130 represents the only solution to Eq. 2-129 
[16,24] . 

Thus the "inverse" filter represents that h r (») function 
which provides the best restoration possible for any com- 
bination of constraints by producing a composite imaging 
system point-spread function having the smallest radius of 
gyration. 



CHAPTER 3 


SYSTEM EQUATION SOLUTION TECHNIQUE FOR 
OPTIMUM RESTORATION FILTER 

3.1 Introduction 

The general system of equations defining the restoration 
filter aperture are examined to establish that a global 
minimum to the minimization problem exists , thus assuring a 
unique solution for a given set of constraints. The second- 
order differential equations of the previous chapter are 
reduced to a system of first-order differential equations to 
facilitate numerical solution, a minimization technique 
based upon the solution of the system of first-order 
differential equations subject to the constraint equations 
is described, and the effect of the initial conditions of 
the system of first-order differential equations upon the 
restoration filter aperture is examined. 
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3.2 Existence of Global Minimum to Restoration Filter 

Equation System 

In the following discussion, direct reference will be 
made to the system of equations which specifies the optimum 
restoration filter for the x-axis component of a separable 
aperture. Eg. 2-9 2, 9 3, 7 6b, and 79b. However, the results are 
completely general and may also be applied to the systems 
of equations defining the y-axis component as well as the 
general restoration filter apertures. 

In order to guarantee that the equation system, Eq. 2-92, 
93, 76b, and 79b, has a unique solution, all the defining 
functionals, that is the functional to be minimized as well 
as the constraint functionals, must be convex [14,15] . The 
property of convexity allows the theory of local extrema for 
general nonlinear functionals to become a global theory. 

Since the augmented functional I , Eq. 2-80, is a linear 
combination of the defining functionals, it must also be 
convex. Thus, there exists only one stationary point of the 
augmented functional and consequently Eq. 2-92 ,93 ,76b ,and 79b 
have a unique solution for a given set of constraints. It 
can be shown [22] that quadratic functionals which are 
squares of norms are convex. Since all the defining functionals 
are of this class, convexity is assured. 



As in the previous section/ for convenience/ specific 
reference will be made to the equation system defining the 
x-axis component of a separable restoration filter aperture. 
However the results of this section are equally applicable to 
the radially symmetric restoration aperture. To include the 
results of the iteration technique described in Section 2.6/ 

Eg, 2-92 will be replaced by Eq. 2-123. In order to take 
advantage of the many numerical techniques available for solving 
first-order differential equations/ it is convenient to 
convert Eq. 2-123 to the normal form which is done in the 
following manner. Eq. 2-123 may be written as 





I 
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and 


C<f x > 


l H bx tf x>i V + A lxV 


(3-4) 


Defining the complex quantities in terms of the 
corresponding real and imaginary components, 

H rx< £ x> " H rxr tf x> + 3 H rxi«x> 


A(£ x> = A r <f x ) + j A i (f J£ ) 


and 


B(f x > = B r i£ x ) + j B.(f x ) 


C(f ) = C (f ) + j C.(f ) 
x' r x' J x x 


(3-5) 

(3-6) 

(3-7) 

(3-8) 


Eq. 3-1 may be rewritten as 

lH- r (f x ) + j H-.(f x >] f [A r (f x )+j A i( f x )] 

• ‘ H ixr< f x> + 3 H ixi( f x>l + ' B r (f x> + 3 W 

• ' H rxr< £ x) + 3 H rxi< £ ,” “ V £ x> + 3 W <3 - 9 > 

or by collecting the real and imaginary components Eq. 3-9 
may be written as 


H iA r < £ x> + V £ x> H ixr !£ x> ‘ V £ x )H ixi (£ x> 


+ B r< £ x 1B rxr< £ x> " B i < £ x> H rxi < £ x> ~ C r (£ x> 

+ 3 IH^ i( f x ) + V £ *»W £ x> 


+ B i< £ x )H rxr l£ x> + B r (£ x )H rxi (£ x> ’ c i (£ x 1! “ °' 


(3-10) 


Eq. 3-10 may be separated into two differential equations 
[16] formed by the real and imaginary components of Eq. 3-10 


+ V £ x> H ixr (f x> - V f x> H ixi tf x> 


- V f x> H r xr [£ x> - V £ x )H rxi‘ £ x> - C r [£ x> = 0 


(3-1 


and 


H”.(f ) + A.(fjH* (f ) + A (f )H' . (f ) 
rxi x' i' x' rxr x' r x rxi x 


+ B. (f )H (f ) + B (f )H . (f ) - C . (f ) = 0. 
i x rxr x r x rxi x i x 


(3-1 


Eq. 3-11 and 12 may be reduced to a system of first- 
order differential equations [16] by introducing the 
variables. 


H. 

(f ) 

= H 

(f ) 

1 

X 

rxr 

X 


(f ) 

= H* 

(f ) 

2 

X 

rxr 

X 

H- 

(f ) 

= H . 

(f ) 

3 

X 

rxi 

X 

H 4 

< £ x> 

s H rxi 

< £ x> 


(3-13) 

By substituting Eq. 3-13 into Eq. 3-11 and 12, the 
following system of first-order differential equations is 
formed , 

H i (f x> - W 


a 2 (f x> “ -V £ x>V f x> - VW f x> 

+ B. (f )H, (f ) + A. (f ) H . (f ) + C <f ) 
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H 3< f x> = V f x ! 

«4 tf x> = - V f x )H l (f x> - WW 

- VW f x> - V f x )H 4< f x> + C i' f x>- < 3 ' 14 > 

where 

H rx (f x ! " H l< f x> + * V f x>- ' 3 - 15 > 

Thus the second-order differential equation represented by 
Eq. 2-123 has been reduced to a system of first-order 
differential equations, Eq. 3-14, 

The technique for determining a solution to the system 
of equations defined by Eq. 3-14, 2-93, 76b , and 79b is 
based upon treating the constraint equations, Eq. 2-93, 76b, 
and 79b, as a system of nonlinear equations and using a 
modified version of the subprogram SECANT [23] which is 
designed to handle a system of nonlinear equations. From 
an initial estimate for the Lagrange multipliers, a 
modified version of the subprogram DHPCG [12] , a subroutine 

designed to solve a system of first-order linear 
differential equations, is used to solve Eq. 3-14 for 
H (f ), The constraint equations, Eq. 2-93, 76b, and 79b, 
are then checked and any resulting error forms the 
basis for another iteration of SECANT and subsequent choice 
of Lagrange multipliers. Although SECANT only converges 
linearly, no derivatives of the functions defining the 
constraint equations are required. 
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Since SECANT has control over DHPCG by varying the 
Lagrange multipliers , it is important to re-examine Eq. 3-14 
to determine the conditions for which no solution to this 
equation exists. It may be shown [16,24] that as long as 
the coefficients of Eq. 3-14, that is A , A., B , B . , C , 
and are continuous functions, then for a given set of 
initial conditions one and only one set of solutions exist. 
Thus theoretically for any set of finite value Lagrange 
multipliers, a solution to Eq, 3-14 should exist. 

To complete the solution for the optimum x-axis 
component of the restoration aperture, appropriate initial 
conditions for Eq. 3-13 must be supplied and are defined 
as. 


H l<°> “ H rxr< 01 


V 0) " H rcr <0> 


h 3 (0) - H^CO) 


H 4 t0 > ‘ 1 W°> 


For h (x ) to be a real function , 

Xa 


cos 


B rxr<V “ J h rx< x > 
H i*r< f x> - 2 " l x \x (x! 


H rxi (f x ) = / h rx (x> sin 

— CO 


then 

2irf x dx 

X 

sin 2irf x dx 

A 

2vf x dx 


(3-16) 


I 
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and 




From Eg. 3-17, 


= 2tt ./ x h (x) cos 2TTf x dx 


x 


(3-17) 


H ixr <0 > = 0 = H rxi< 01 - 


(3-18) 


The remaining initial conditions may he used to 

control the shape of (x) . For h „ (x) to be an even 
*■ rx rx 

function about the x-axis origin, then 


H rxr' 01 * 0 


and 




(3-19) 


To preserve the DC gain, or average value, of the blurring 
system. 


H_(0) = 


(3-20) 


rxr 


Ret H bx (0) ] 


where Re[.] denotes the real part of the complex argument. 
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CHAPTER 4 

EXPERIMENTAL VERIFICATION OF FUNDAMENTAL PROPERTIES OF 
OPTIMUM RESTORATION FILTER 

4 . 1 Introduction 

This chapter presents experimental results demonstrating 
the restoration filter performance as a function of the con- 
straint parameters and from its application to a series of 
blurred test patterns. These results are based upon the 
assumption of a Gaussian-shaped blurring aperture. Since 
the principle application for this restoration technique 
was the resolution enhancement of the multispectral scanner 
system used for ERTS data collection, and since no specific 
information concerning the shape of the multispectral scanner 
aperture was available, a Gaussian model was chosen. This 
assumption was based upon two primary considerations. 

Because of the advantages in correcting a separable aperture 
as previously described in section 2.**, the first consid- 
eration for choosing a Gaussian blurring aperture was that 
such a model, if chosen to be symmetric about each axis, 
could also provide a separable blurring function which was 
also radially symmetric, thus approximating the radial 
symmetry of the actual scanner aperture. The second con- 
sideration for choosing a Gaussian blurring aperture model 
was for mathematical convenience. Since the frequency 
spectrum of a Gaussian blurring aperture is also Gaussian, 
this spectrum may be explicitly computed, thus eliminating 
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any errors which would have otherwise been introduced into 
the computation of the blurring aperture spectrum, or 

as defined in Eq. 2-92 or 2-96 respectively, by a sub- 
routine such as FORT. 

Based upon a Gaussian blurring aperture model, an 
investigation of the performance and shape of the restoration 
filter as a function of its fundamental parameters and con- 
straints is made. Also the composite system performance as 
a function of blurring aperture error is examined. 

In addition, the blurring effects of a specific Gaussian 
aperture on a series of test patterns are examined to 
qualitatively demonstrate the degree of resolution enchance- 
ment possible for a given set of constraints. Also the trun- 
cation error effects for the "classical" inverse and con- 
strained restoration filters are compared. 

4.2 Parameter and Constraint Control of Restoration Filter 

Performance 

The following is an investigation of the controlling 
effect demonstrated by the defined parameters and constraints 
upon the restoration filter's shape and performance. The 
constraint variations considered were directed by two 
primary goals: the first to investigate how parameter 

variations affect the restoration filter, and the second 
to apply the restoration filter to a specific enhancement 
problem which will be discussed in the next chapter. The 
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restoration filters for three Gaussian blurring apertures 
having radii of gyration of 1, 3, and 5, respectively, will 
be examined and compared to the unconstrained or classical 
"inverse" restoration filter. Since a separable blurring 
aperture is assumed, the following results may be interpreted 
as either x-axis solutions to the system of operations 
defined by Eg. 2-123, 121, 93, 76b, and 79b, or as y-axis 
solutions defined by a similar equation system where x- 
dependent parameters are replaced by corresponding y- 
dependent parameters , 

Figure 4.1 shows a Gaussian shaped blurring aperture 
with a radius of gyration of one. The frequency spectrum 
of the classical "inverse" restoration filter, is shown in 
Figure 4.2. It should be noted that this spectral data was 
obtained not by inverting the spectrum of the blurring 
aperture of Figure 4.1 but rather by imposing the conditions 
stated in section 2.7 which define the classical "inverse" 
filter in terms of the parameters of the constrained resto- 
ration filter, Eq. 2-123, 121, 93, 76b, and 79b, and was 
performed as a test to verify the accuracy of the overall 
numerical algorithm described in Chapter 3. The correspond- 
ing point-spread function of the classical "inverse" 
restoration filter is shown in Figure 4.3. As a further 
check on the overall accuracy of the algorithm, the composite 
system point-spread function, representing the convolution 
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Figure 4,1 Gaussian Blurring Aperture 
Having a Radius of Gyration of One 
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Figure 4.2 Classical "Inverse" Restoration 
Filter Frequency Spectrum for Blurring Aperture 
of Figure 4.1 
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Figure 4.3 Classical "Inverse" Restoration 
Filter Point-Spread Function for Blurring Aperture 
of Figure 4.1 



Figure 4.4 Composite System Point- Spread Function 
Resulting from Correcting the Blurring Aperture of 
Figure 4.1 with the Restoration Function of Figure 
4.3 


of the blurring aperture of Figure 4.1 with the "inverse" 
restoration function of Figure 4.3 is shown in Figure 4.4. 
As previously stated in section 2.2, the classical "inverse 
filter will produce an impulse for the composite system 
point- spread function. Thus Figure 4.4 should represent a. 
digital approximation to an impulse. Since the base line 
of the function shown in this figure should be zero, the 
difference, approximately 7 x 10 ^ , represents a measure 
of the round-off error within the algorithm. 

Figure 4.5 represents the restoration filter frequency 
spectrum for the blurring aperture of Figure 4.1 when the 
constraint, K^, is employed to limit the radius of gyration 
of the resulting restoration filter point-spread function 
to 0.S that of the unconstrained or classical "inverse" 
filter. Figure 4.6 shows the resulting restoration filter 
point-spread function having a radius of gyration of 0.5 
that of the classical "inverse" restoration function shown 
in Figure 4.3. The resulting composite system point-spread 
function produced by convolving the blurring aperture of 
Figure 4.1 with the restoration point-spread function of 
Figure 4.6, is shown in Figure 4.7. The resulting radius 
of gyration of the composite system point-spread function 
is approximately 67% of that of the blurring aperture. 
Particular attention, however, should be paid to the signif 
icant secondary responses which are approximately 10% of 
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Figure 4.5 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.1 for K = 0.5 
Unconstrained Value with No Criterion Function 
Iteration 
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Figure 4.6 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.5 
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Figure 4.8 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.1 for = 0.5 
Unconstrained Value and Nine Iterations with an 
Exponential Criterion Function 
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the magnitude of the primary lobe. As described in section 
2.6/ these secondary responses are highly undesirable since 
ghost images in the corrected image may result. 

By applying the iteration technique described in section 

2.6 to the restoration function previously described by 
Figure 4. 5- 4. 7 and using nine iterations with an exponential 
criterion function, the resulting restoration filter fre- 
quency spectrum for the same constraint on is shown in 
Figure 4.8. The corresponding restoration filter point- 
spread function is shown in Figure 4.9 and has the same 
radius of gyration as the restoration function shown in 
Figure 4.6. The composite system point-spread function/ 
resulting from convolving the function of Figure 4.1 with 
that of Figure 4.9, is shown in Figure 4.10. By comparing 
the composite system point-spread function shown in Figure 

4.7 with that of Figure 4.10, the effectiveness of the 
iteration technique to suppress secondary side-lobe re- 
sponses is evident. The radius of gyration of the composite 
system point-spread function of Figure 4.10 is approximately 
74% of that the blurring aperture. Thus, as expected, the 
restoration function resulting from the iteration technique 
is somewhat less effective for image enhancement. It 
should be noted that although nine iterations were used, 
there is little change in the restoration function after 
about three iterations and that the iteration technique 
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Figure 4.9 Restoration Filter Point- Spread 
Function Having Spectrum of Figure 4.9 
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Figure 4.10 Composite System Point-Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.1 with the Restoration 
Function of Figure 4.9 
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appears to converge rapidly. 

Figure 4.11 through 4.16 and 4.17 through 4.22 represent 
the resulting restoration and composite system functions 
when the constraint in the radius of gyration of the resto- 
ration filter, is further reduced to 0.2 and 0.1 of the 
radius of gyration of the "inverse" filter, respectively. 

The radii of gyration of th^ composite system point-spread 
functions shown in Figure 4.13, 4.16, 4.19, and 4.22 are 81, 
97, 89 and 90 percent of the radius of gyration of the 
blurring aperture of Figure 4.1, respectively. 

The preceding figures demonstrated the effect of the 
constraint on the radius of gyration of the restoration 
filter in a noiseless environment for a Gaussian blurring 
aperture with a radius of gyration of one. To demonstrate 

the influence of noise on the restoration filter, an 

-5 

additive white noise with a spectral density of 10 
2 

volts /Hz was assumed. In addition, the output noise power 
constraint, K^, was chosen to equal approximately 0.149% of 
the corresponding output noise power of the "inverse" filter. 
This particular value was chosen to approximate a signal- 
to-noise ratio of 20 dB when applied to the multispectral 
scanner data to be described in a later section. Figure 
4.23 and 4.24 show the frequency spectra of the restoration 
filters when the constraint on the radius of gyration of 
the restoration filter point-spread function equals 0.5 
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Figure 4.11 Restoration Filtex Frequency Spectrum 
for Blurring Aperture of Figure 4.1 for K 1 = 0.2 
Unconstrained Value with Wo Criterion Function 
Iteration 



Figure 4.12 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.11 


] 




im 


NORMALIZED DISTANCE 

Figure 4.13 Composite System Point-Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.1 with the Restoration 
Function of Figure 4.12 
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Figure 4.14 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.1 for K, = 0.2 
Unconstrained Value with Ten Iterations Using an 
Exponential Criterion Function 
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Figure 4.15 Restoration Filter Point- Spread 
Function Having Spectrum of Figure 4.14 



Figure 4.16 Composite System Point- Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.1 with the Restoration 
Function of Figure 4.15 
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Figure 4.17 Restoration Filter Frequency 
Spectrum for Blurring Aperture of Figure 4.1 
for K.*0.1 Unconstrained Value With No Criterion 
Function Iterations 
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Figure 4.18 Restoration Filter Point- Spread 
Function Having Spectrum of Figure 4.17 
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Figure 4.19 Composite System Point-Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.1 with the Restoration 
Function of Figure 4.18 
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Figure 4.20 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.1 for K, =0.1 
Unconstrained Value with Ten Iterations Using an 
Exponential Criterion Function 
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Figure 4.21 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.20 
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Figure 4.22 Composite System Point-Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.1 with the Restoration 
Function of Figure 4.21 
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Figure 4.23 Restoration Filter Spectrum for 
Blurring Aperture of Figure 4.1 for K = 0.5 and 
K „ = 0.00149 (S/N » 20dB) Unconstrained Values__ 
for a White Noise Power Spectral Density of 10** 5 
volts 2 /Hz 
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Figure 4.25 Gaussian Blurring Aperture Having 
a Radius of Gyration of Three 





80 


and 0.2 times the value of the "inverse" filter, respective- 
ly. 

The following discussion will be with reference to a 
Gaussian blurring aperture with a radius of gyration of 
three, as shown in Figure 4.25. The restoration filter 
frequency spectrum, point-spread function, and the composite 
system point-spread function for the classical "inverse" 
restoration filter are shown in Figure 4.26, 4.27 and 4.28, 
respectively. When the radius of gyration of the resto- 
ration filter point-spread function is constrained to 0.5 
of the unconstrained value, or value of the "inverse" filter, 
the resulting restoration and composite system functions 
with and without criterion function iterations are shown 
in Figure 4.29 to 4.34. Figure 4.35 to 4 U 40 show correspon- 
ding results when the radius of gyration of the restoration 
filter is further constrained to 0.2 of the value of the 
"inverse" filters. The radii of gyration of the composite 
system point-spread functions of Figure 4.31, 4.34, 4.37 and 
4.40 are 14.3, 19.8, 88.5, and 89.3 percent of the radius 
of gyration of the blurring aperture of Figure 4.25. 

Figure 4.41-4.52 demonstrate the effect of a white 
noise power spectral density upon the restoration function 
for the same input and output signal-to-noise ratios as 
employed with the filter previously considered. Figure 
4.41-4.46 . show the frequency spectrum of the restoration 
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Figure 4.26 Classical "Inverse" Restoration 
Filter Frequency Spectrum for Blurring Aperture 
of Figure 4.25 
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Figure 4.27 Classical "Inverse" Restoration 
Filter Point- Spread Function for Blurring Aperture 
of Figure 4,25 
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Figure 4.28 Composite System Point- Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.25 with the Restoration 
Function of Figure 4.27 
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Figure 4.29 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.25 for K. = 0.5 
Unconstrained Value with No Criterion Function 
Iteration 
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Figure 4.30 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.29 
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Figure 4.31 Composite System Point-Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.25 with the Restoration 
Function of Figure 4.30 
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Figure 4.32 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.25 for K, = 0.5 
Unconstrained Value and Four Iterations with ai? 
Exponential Criterion Function 
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Figure 4.33 "Restoration Filter Point- Spread 
Function Having Spectrum of Figure 4.32 
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Figure 4.35 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.25 for K_ = 0.2 
Unconstrained Value with No Criterion Function 
Iteration 
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Figure 4.36 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.35 
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Figure 4.37 Composite System Point- Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.25 with the Restoration 
Function of Figure 4.36 
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Figure 4.38 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.25 for K. = 0.2 
Unconstrained Value and Four Iterations with an 
Exponential Criterion Function 
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Figure 4.39 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.38 
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Figure 4.40 Composite System Point-Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.25 with the Restoration 
Function of Figure 4.39 
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Figure 4.41 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.25 for K 1 =0.5 and 
K = 1.61 x 10 -7 (S/N ~ 20dB) Unconstrained Values 
for a White Noise Power Spectral Density of 10*5 
volts^/Hz and Four Iterations Using an Exponential 
Criterion Function 
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Figure 4.42 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.41 



Figure 4.43 Composite System Point- Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.25 with the Restoration 
Function of Figure 4.42 
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Figure 4.44 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.25 for K. = 0.2 
and K„ = 1.61 x 10“ ^ (S/N « 20dB) Unconstrained 
Values for a White Noise Power Spectral Density of 
10-5 volts2/Hz and Four Iterations Using an Exponential 
Criterion Function 
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Figure 4.45 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.44 


NORMALIZED DISTANCE 

Figure 4.46 Composite System Point-Spread 

Function Resulting from Correcting the Blurring 

Aperture of Figure 4.25 with the Restoration 

Function of Figure 4.45 


0.00 16.00 3e.00 46.00 

NORMALIZED FREQUENCY 

Figure 4.47 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.25 for K, = 0.5 and 
K = 1.61 x 10"" 8 (S/N s 20dB) Unconstrained Values 
for a White Noise Power Spectral Density of 10~ 4 
volts 2 /Hz and Four Iterations Using an Exponential 
Criterion Function 
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filter and the point-spread functions for both the resto- 
ration filter and the composite system for both = 0.5 

and 0.2 of the unconstrained value for a white noise power 

-5 2 

spectral density of 10 volts /Hz. Figure 4.47-4.52 show 

the effect of increasing the noise power spectral density 
-4 -3 2 

to 10 and 10 volts /Hz for = 0.5 of the unconstrained 

value. The radii of gyration of the composite system point- 
spread functions shown in Figure 4.43, 4.46, 4.49 and 4.52 
are 60.2, 88.5, 66.2, and 76.9 percent of the radius of 
gyration of the blurring aperture of Figure 4.25. 

Thr point-spread function for a Gaussian blurring 
aperture having a radius of gyration of five is shown in 
Figure 4.53. Figure 4.54 shows the resulting "inverse” 
restoration filter frequency spectrum, while the correspon- 
ding "inverse" restoration filter point-spread function is 
shown in Figure 4.55. The overall composite system point- 
spread function, resulting from correcting the blurring 
aperture of Figure 4.53 with the "inverse" restoration 
function of Figure 4.55, is shown in Figure 4.56. 

When the radius of gyration of tne point-spread function 
of the restoration filter is constrained to 0.5 of that of 
the "inverse" filter, the resulting restoration filter 
frequency spectrum, point-spread and composite system point- 
spread function are shown in Figure 4.57-4.59, respectively. 
By further reducing this constraint to 0.2 of that of the 
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Figure 4.48 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.47 



Figure 4.49 Composite System Point-Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.25 with the Restoration 
Function of Figure 4.48 
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Figure 4.50 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.25 for K. = 0.5 and 
K„ = 1.61 x 10~9 (S/N s 20dB) Unconstrained Values 
for a White Noise Power Spectral Density of 10“^ 
volts 2 /Hz and Two Iterations Using an Exponential 
Criterion Function 
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Figure 4.51 Restoration Filter Point— Spread 
Function Having Spectrum of Figure 4.50 




NORMALIZED DISTANCE 

Figure 4.52 Composite System Point-Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.25 with the Restoration 
Function of Figure 4.51 
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Figure 4.53 Gaussian Blurring Aperture Having 
a Radius of Gyration of Five 
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Figure 4.54 Classical "Inverse” Restoration 
Filter Frequency Spectrum for Blurring Aperture 
of Figure 4.53 
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Figure 4.55 Classical "Inverse" Restoration 
Filter Point- Spread Function for Blurring Aperture 
of Figure 4.53 
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Figure 4.56 Composite EJystem Point-Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.53 with the Restoration 
Function of Figure 4.55 
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Figure 4.57 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.53 for K. = 0.5 
Unconstrained Value and Two Iterations with an 
Exponential Criterion Function 
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Figure 4.58 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.57 


1.000-1 


.600 ■ 


0 .000 ■ 


g(x) 


— j 1 1 1 

0.Q0 16.00 32.00 40.00 64.00 

NORMALIZED DISTANCE 

Figure 4.59 Composite System Point-Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.53 with the Restoration 
Function of Figure 4.58 
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"inverse" filter, the resulting restoration filter frequency 
spectrum, point- spread function and composite system point- 
spread function are shown in Figure 4.60-4.62, respectively. 

Figure 4.63-4.65 show the restoration filter frequency 
spectrum, point-spread function, and composite system point- 
spread function when the constraint on the radius of gyration 
of the restoration filter point-spread function is 0.1 of 
that of the "inverse" filter. The radii of gyration of the 
composite system point-spread functions shown in Figure 4.59, 
4.62, and 4.65 are 15, 21, and 50 percent of the radius of 
gyration of the blurring aperture of Figure 4.53. 

The preceding results corresponding to the blurring 
aperture of Figure 4.53 were for a noiseless environment. 

The following results are based upon two types of noise 
power spectral densities: white and exponentially increasing 

functions. The white noise power spectral density will be 
considered first. Results for four spectral density levels 
are shown: 10 ^ , 10~^, 5x10 and 10 2 volts 2 /Hz. The 

output noise level constraint in the restored image is 
adjusted for a signal- to-noise ratio of approximately 20 dB. 
Figure 4.66-4.68 show the restoration filter frequency 
spectrum, point-spread function, and composite system 
point-spread function for K^, the restoration filter point- 
spread function radius of gyration constraint equal to 0.5 
that of the "inverse" filter, and I< 2 , the output noise 
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Figure 4.60 Restoration Filter Frequency 
Spectrum for Blurring Aperture of Figure 4.53 
for K,=0.2 Unconstrained Value and Two 
Iterations With An Exponential Criterion 
Function 
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Figure 4.61 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.60 
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Figure 4.62 Composite System Point-Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.53 with the Restoration 
Function of Figure 4.61 
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Figure 4.63 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.53 for K. = 0.1 
Unconstrained Value and Ten Iterations witn an 
Exponential Criterion Function 
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Figure 4.64 Restoration Filter Point- Spread 
Function Having Spectrum of Figure 4,63 
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Figure 4.66 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.53 for K 1 = 0.5 and 
K„ — 2.10 x 10“ i 3 (S/N ~ 20dB) Unconstrained Values 
for a White Noise Power Spectral Density of 10“ 5 
volts 2 /Hz and Ten Iterations with an Exponential 
Criterion Function 
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Figure 4.67 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.66 



Figure 4.68 Composite System Point-Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.53 with the Restoration 
Function of Figure 4.67 
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power constraint, adjusted for a 20 dB output signal-to- 
noise ratio when the noise power spectral density level is 
10” 5 volts ^/Hz. Figure 4.69-4.71 show similar results when 
= 0.2 that of the "inverse” filter. Figure 4.72-4.74, 
4.75-4.77, and 4.78-4.80 show the restoration filter fre- 
quency spectra point-spread functions, and composite system 

point-spread functions for = 0.5 that of the "inverse” 

-4 

filter and for noise power spectral densities of 10 , 

5xl0~ 4 , and 10 -3 volts 2 /Hz, respectively. 

The exponentially increasing noise power spectral 
densities are shown in Figure 4.81, 4.85, and 4.89. Figure 
4.82-4.84, 4.86-4.88, and 4.90-4.92 show the corresponding 
restoration filter frequency spectra, point-spread functions 
and composite system point-spread functions, respectively 
for = 0.5 that of the "inverse" filter and K 2 adjusted 
for a signal-to-noise ratio of 20 dB, for each noise 
spectrum. 

Tables 4.1 and 4.2 summarize the basic restoration filter 
and composite system performance parameters for the con- 
straint and noise power spectral density variations described 


in this section. 
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Figure 4.69 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.53 for K. = 0.2 and 
= 2.10 x 10~-L3 (S/N ~ 20dB) Unconstrained Values 
for a White Noise Power Spectral Density of 10”® 
volts 2 /Hz and Four Iterations with an Exponential 
Criterion Function 
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Figure 4.70 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.69 
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Figure 4.71 Composite System Point Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.53 with the Restoration 
Function of Figure 4.70 



Figure 4.72 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.53 for K. = 0.5 and 
K_ = 2.10 x 10“-*-4 (S/N 3 20dB) Unconstrained Values 
for a White Noise Power Spectral Density of 10“ 4 
volts 2 /Hz and Four Iterations with an Exponential 
Criterion Function 
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Figure 4.73 Restoration Filter Point- Spread 
Function Having Spectrum of Figure 4.72 
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Figure 4.74 Composite System Point-Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.53 with the Restoration 
Function of Figure 4.73 
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Figure 4.75 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.53 for K.. = 0.5 and 
K_ = 1.05 x 10”15 (s/N « 20dB) Unconstrained Values 
for a White Noise Power Spectral Density of 5 x 10~4 
volts^/Hz and Four Iterations with an Exponential 
Criterion Function 



Figure 4.76 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.75 
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Figure 4.77 Composite System Point- Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.53 with the Restoration 
Function of Figure 4.76 



Figure 4.78 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.53 for K. = 0.5 and 
K = 2.10 X 10” 15 (S/M ~ 20 dB) Unconstrained Values 
for a White Noise Power Spectral Density of 10“ 3 
volts 2 /Hz and Four Iterations with an Exponential 
Criterion Function 
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Figure 4.79 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.78 
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Figure 4.80 Composite System Point- Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.53 with the Restoration 
Function of Figure 4.79 
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Figure 4.81 N 
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Figure 4.82 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.53 for Noise Power 
Spectral Density of Figure 4.81 and K_ = 0.5 Uncon~ 
strained Value and K„ = 5 IS/N z 20dBj for Four 
Iterations with an Exponential Criterion Function 
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Figure 4.83 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.82 


0.00 

NORMALIZED DISTANCE 

Figure 4.84 Composite System Point- Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.53 with the Restoration 
Function of Figure 4.83 
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Figure 4,85 Noise Power Spectral Density for 
$ (f) = 10”5 £ v^/Hz 
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Figure 4.86 Restoration Filter Frequency Spectrum 
for Blurring Aperture of Figure 4.53 for Noise 
Power Spectral Density of Figure 4.85 and K, = 0.5 
Unconstrained Value and K 2 = 5 (S/N ~ 20dB) x for Four 
Iterations with an Exponential Criterion Function 
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Figure 4.88 Composite System Point-Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.53 with the Restoration 
Function of Figure 4.87 
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Figure 4.91 Restoration Filter Point-Spread 
Function Having Spectrum of Figure 4.90 
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Figure 4.92 Composite System Point-Spread 
Function Resulting from Correcting the Blurring 
Aperture of Figure 4.53 with the Restoration 
Function of Figure 4.91 


Table 4.1: Basic Restoration Filter and Cd 


blurring aperture radius of gyration 


K 1 constraint as a fraction of the unconstrained value 


K 2 constraint for approximate 20dB output signal-to-noise ratio 


input signal-to-noise ratio for specified white noise power 
spectral density amplitude 


criterion function iteration 


composite system radius of gyration as percentage of blurring 
aperture radius of gyration 


maximum peak of secondary oscillations as percentage of primary 
lobe in composite system point-spread function 


blurring aperture radius of gyration 

constraint as a fraction of the unconstrained value 

x 2 constraint for approximate 20dB output signal-to-noise ratio 


input signal-to-noise ratio for specified white noise power 
spectral density amplitude 


criterion function iteration 


composite system radius of gyration as percentage of blurring 
aperture radius of gyration 


maximum peak of secondary oscillations as percentage of primary 
lobe in composite system point-spread function 


blurring aperture radius of gyration 


constraint as a fraction of the unconstrained value 


k 2 constraint for approximate 20dB output signal-to-noise ratio 


input signal-to-noise ratio for specified white noise power 
spectral density amplitude 


criterion function iteration 


composite system radius of gyration as percentage of blurring 
aperture ra^us of gyration 


maximum peak of secondary oscillations as percentage of primary 
lobe in composite system point-spread function 
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Table 4.1/ cont 


blurring aperture radius of gyration 
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K 2 constraint for approximate 20dB output signal- to-noise ratio 
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maximum peak of secondary oscillations as percentage of primary 
lobe in composite system point- spread function 
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Table 4.2: Basic Restoration Filter and 

Performance Parameters As A 
and Exponential Noise Spectz 


blurring aperture radius of gyration 
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lobe in composite system point-spread function 


0.5 


5 

4 <£) - 10“ 6 f 1-8S v 2 /Hz 

nn 

29.57dB 

none 4 r exponential 


60 


66 


-1.38 


-0.30 


ORIGINAL P^EIS 

Op POOR QtTALTU 


HOLD OUT FRA 11^ 


n Filter and Composite System 
meters As A Function of Constraints 
floise Spectral Density Variations 


5 



0*5 


0.5 


5 


5 

v 2 /llz 

<5 (f) = 10“ 5 f V 2 /Hz 

nn 

nn 

= l(f 5 f 1 - 85 V 2 /HZ 


47* 72dB 


19.57dB 

exponential 

none 4, exponential 

none 

4, exponential 

66 

SZ 73 

89 

94 

"0*36 

-1,42 -0.39 

-0,62 

-3,4 x 10 -3 




*" 


FQLD OlTg ER AME} ^ 


121 


4.3 The Effect Of Error Between Actual And Assumed 

Gaussian Blurring Apertures' On Composite System 

Performance 

The results of the preceding section assumed that the 
blurring aperture function was explicitly known. When this 
function is not known, the composite system performance 
based upon an erroneous blurring function can deviate 
significantly from the optimal response. The purpose of 
this section is to demonstrate the effect upon the composite 
system point-spread function produced when a restoration 
filter computed for a Gaussian blurring aperture of spe- 
cified radius of gyration is used to correct data blurred 
by a Gaussian aperture having a different radius of gyration. 

Figure 4.93-4.96 show both the composite system fre- 
quency spectra and point-spread functions when restoration 
filters computed for Gaussian blurring apertures having 
radii of gyration of three and five, respectively, are used 
to correct a Gaussian blurring aperture having a radius of 
gyration of one. Figure 4.97 and 4.98 show the composite 
system frequency sepctrum and point-spread functions when 
a restoration filter computed for a Gaussian blurring aper- 
ture having a radius of gyration of five is used to correct 
a Gaussian blurring aperture having a radius of gyration 
of three. These preceding mismatch combinations of actual 
and assumed blurring apertures correspond to an overcorrec- 
tion of the data as demonstrated by the peaks in the com- 
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Figure 4.93 Composite System Frequency Spectrum 
Resulting from Correcting the Blurring Aperture of 
Figure 4.1 with the Restoration Filter of Figure 
4.45 


NORMALIZED DISTANCE 

Figure 4.94 Composite System Point- Spread 
Function Having Spectrum of Figure 4.93 
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Figure 4.95 Composite System Frequency Spectrum 
Resulting from Correcting the Blurring Aperture of 
Figure 4.1 with the Restoration Filter of Figure 
4.7 
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Figure 4.96 Composite System Point- Spread 
Function Having Spectrum of Figure 4.95 
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Figure 4.27 Composite System Frequency Spectrum 
Resulting from Correcting the Blurring Aperture of 
Figure 4.25 with the Restoration Filter of Figure 
4.70 



Figure 4.98 Composite System Point- Spread 
Function Having Spectrum of Figure 4.97 
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posite system frequency spectra of Figure 4.93, 4.95, and 
4.97, The resulting composite system point-spread functions, 
as shown in Figure 4.94, 4.96, and 4.98 are highly oscilla- 
tory - the greater the assumed blurring aperture relative 
to the actual, the greater the oscillations in the composite 
system point- spread functions. 

Figure 4.99-4.104 show corresponding composite system 
spectra and point-spread functions for the opposite situation, 
where the radius of gyration of the actual blurring aperture 
is larger than that of the assumed aperture. These blurring 
aperture mismatch combinations correspond to undercorrection 
of the data. 

Of the two types of mismatch conditions considered, the 
least desirable is the overcorrection condition; since it 
would introduce ghosts into the "corrected” image. In con- 
trast, the undercorrected condition produces a stable system 
response which would not introduce any oscillatory error; 
however, less than optimal image enhancement could be 
expected. 

The recognition of the overcorrection condition by a 
peak in the frequency spectrum of the "corrected" data or 
by the presence of oscillations about transients in the 
"corrected" image could provide a criterion for determining 
when a blurring aperture match has been achieved for data 
where the radius of gyration of the actual Gaussian- shaped 
aperture is unknown. 
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Figure 4.99 Composite System Frequency Spectrum 
Resulting from Correcting the Blurring Aperture of 
Figure 4.25 with the Restoration Filter of Figure 
4.24 
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Figure 4.100 Composite System Point- Spread 
Function Having Spectrum of Figure 4.99 
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Figure 4.101 Composite System Frequency Spectrum 
Resulting from Correcting the Blurring Aperture of 
Figure 4.53 with the Restoration Filter of Figure 
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Figure 4.103 Composite System Frequency Spectrum 
Resulting from Correcting the Blurring Aperture of 
Figure 4.53 with the Restoration Filter of Figure 
4.45 
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Figure 4.104 Composite System Point-Spread 
Function Having System of Figure 4.103 
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4.4 Application of the Restoration Filter to linages Blurred 

by a Known Gaussian Aperture 

The data presented in this section was chosen to qual- 
itatively demonstrate the fundamental properties of the 
optimum restoration function as previously described analyt- 
ically in Chapter 2 and experimentally presented in section 2 
of this chapter. Figure 4-105 thru 107 demonstrate the 
resolution improvement possible when a given test pattern 
initially blurred by a known aperture function is reconstructed 
under essentially noiseless conditions. It should be men- 
tioned that round-off errors within the restoration program, 
resulting principally from FORT, are primarily responsible 
for the noise present. The input image to the blurring 
system is shown in Figure 4-105. This pattern is described 
by a two dimensional array of 512 x 512, as are all the 
images shown and is composed of a series of rectangular 
wedges whose periodicity ranges from 16 pixels/cycle to 2 
pixels/cycle in 2 integer increments. This frequency 
variation was chosen to provide a simple qualitative com- 
parison of the resolution of the blurred image to that of 
the restored image. In addition, the 8 increment grey level 
pattern was included to demonstrate any changes to the 
dynamic range of the blurred and restored images. Figure 4- 
106 represents the output of an imaging system having a 
Gaussian blurring aperture with a radius of gyration of 3, 
as shown in Figure 4-25. It is evident from this figure 
that the constrast of even the lowest frequency wedge has 



Figure 4.105 Resolution Bar Test Pattern 



Figure 4.106 Figure 4.105 Blurred by a Gaussian 
Aperture with a Radius of Gyration of Three 




been significantly reduced and that the highest frequency 
present in this pattern has a period of approximately 12 
pixels/cycle. By applying the restoration function of 
Figure 4-42 to the data of Figure 4-106 a significant im- 
provement in picture quality is produced as shown, in Figure 
4-107. The highest frequency present in the restored image 
has a period between 6 and 4 pixels/cycle, thus providing 
more than a 2 to 1 resolution improvement. 

The restoration function chosen for this resolution 
improvement test provided a high output signal-to-noise 
ratio, over 40dB, for the round-off error noise introduced 
by the restoration program. Had the round-off error noise 
been less, a more powerful restoration filter could have 
been employed giving a greater improvement in resolution for 
the same output signal-to-noise ratio. It should also be 
noted that had no round-off error noise been present, the 
classical "inverse” filter would have provided perfect 
restoration of the blurred image. 

The test results shown in Figure 4-108 thru 111 demon- 
strate the "ghost" effect produced by restoring a blurred 
image with a restoration function not specifically designed 
to control secondary oscillations about the origin of the 
composite system point- spread function. This control, as 
described in Chapter 2, was obtained by application of the 
criterion function iteration procedure. The input image. 
Figure 4-108, represents a bar pattern having a period of 









Figure 4.107 Figure 4.106 Restored by the 
Restoration Function of Figure 4.42 



Figure 4.108 Low Frequency Square Bar Pattern 
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100 pixels/cycle. Figure 4-109 represents the output of an 
imaging system having a Gaussian blurring aperture with a 
radius of gyration of 3, the same as used in the previous 
set of figures. Figure 4-110 represents a restoration of 
Figure 4-109 with the noniterated version of the restoration 
function of Figure 4-42. A "ghost" image or oscillation on 
each side of the discontinuities cf the bar pattern is 
apparent. However , when the iterated restoration function 
of Figure 4-42 is used to correct the image of Figure 4-109 , 
the resulting image. Figure 4-111, shows no evidence of any 
"ghosts". 

The final fundamental property of the optimum restoration 
filter, and possibly the most significant is the reduction of 
truncation error as compared to that produced by the classi- 
cal "inverse" filter and which is produced by restoration of 
a truncated version of the original blurred image. This 
effect will first be demonstrated by a one-dimensional 
example. Figure 4-112 represents a rectangular test function 
which when blurred by the Gaussian aperture of Figure 4.53 
having a radius of gyration of 5 produces the data of Figure 
4-113. Figure 4-114 is obtained by truncating, or setting 
to zero, the first 10 points of Figure 4-113. By applying 
the classical "inverse" restoration filter of Figure 4-55 
to the truncated data set of Figure 4-114 , the data of Figure 
4-115 results. At no point in the data record does the 
truncation error become small enough to resolve the ' riginal 


I 


I 


I 


I 



Figure 4.109 Figure 4.108 Blurred by a Gaussian 
Aperture with a Radius of Gyration of Three 



Figure 4.110 Figure 4.109 Restored by the Uniterated 
Version of the Restoration Function of Figure 4.42 
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Figure 4.115 "Restored" Function Obtained By 
Correcting the Truncated Data of Figure 4.114 With 
the Classical "Inverse" Restoration Function of 
Figure 4.55 



Figure 4.116 Correction of Truncated Data of 
Figure 4.114 by the Restoration Function of 
Figure 4.67 
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teat function. When the restoration function of Figure 4-67 
is applied to the same truncated data set of Figure 4-114 the 
result, shown in Figure 4-116, is obtained and indicates 
significant error at the edges of the restoration. However, 
because the radius of gyration of the restoration filter 
point-spread function was constrained negligible error is 
present at the center of the data record. 

Figure 4-117 thru 4-120 demonstrate this same effect 
for a two-dimensional example. Figure 4-117 represents a 
bar pattern with a spatial period of 10 pixels/cycle. When 
the Gaussian aperture of Figure 4-25 is used to blur the 
image of Figure 4-117, the result is shown in Figure 4-118. 

By truncating, i.e. setting to zero, all but the center 
256 x 256 portion of Figure 4-118 and applying the classical 
"inverse" restoration filter to the result, the image of 
Figure 4-119 is produced. Although the truncation error is 
smallest at the center of this image, no significant res- 
toration of the original tar pattern is apparent. However, 
when the restoration function of Figure 4-42 is applied to 
the same truncated version of Figure 4-118, the result is 
shown in Figure 4-120. It is clearly evident that the 
truncation error rapidly becomes negligible inside the inter 
256 x 256 portion of this figure and that significant res- 
toration of the original bar pattern is apparent. 


Figure 4.117 High Frequency Square Bar Pattern 



Figure 4.118 Figure 4.117 Blurred by a Gaussian 
Aperture with a Radius of Gyration of Three 
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Figure 4.119 Restoration of a 256 Center Portion 
of Figure 4.118 by the Classical "Inverse" Restoration 
Function of Figure 4.27 



2 

Figure 4.120 Restoration of a 256 Center Portion 
of Figure 4.118 by the Restoration Function of 
Figure 4.42 
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CHAPTER 5 

APPLICATION OF OPTIMUM IMAGE RESTORATION FILTER 
TO ERTS MULTI SPECTRAL SCANNER DATA 

5.1 Introduction 

The primary motivation for this research was the devel- 
opment of an optimal resolution enhancement technique for 
application to remotely sensed multispectral scanner data. 
The fundamental properties of the resulting optimal resto- 
ration filter examined in the previous chapter indicated 
that the shortcomings of the classical "inverse" filter , 
which formed the basis for defining most of these properties 
as discussed in Chapter 2, had been resolved. However, 
until this restoration technique was applied to actual ERTS 
data, it would not be known whether there existed any 
significant deficiencies which had not been previously 
anticipated. 

This chapter describes the application of this resto- 
ration technique to a specific ERTS data set. The problems 
encountered in estimating the multispectral scanner blurring 
aperture, defining the noise present in the data, and per- 
forming the correction to the data are examined. Finally 
conclusions are drawn concerning the quality of the resto- 
ration achieved and the effect of numerical error both in 
the calculation and implementation of the restoration 
function. An IBM 360 model 67 computer was used to compute 
the restoration filters of Chapter 4 and to apply these 


3.42 

corrections to the selected ERTS data set. 

5.2 Estimation of Blurring Aperture 

As mentioned earlier, since no specific information 
concerning the blurring aperture of the multispectral 
scanner was available, a mathematical modal had to be defined. 
As discussed in the introduction to Chapter 4, a Gaussian 
model was ultimately chosen in order to make possible certain 
simplifying assumptions. These assumptions include defining 
a separable blurring process to substantially reduce proces- 
sing cost and time and to permit explicit computation of the 
blurring aperture spectrum, so that numerical round-off 
errors may be reduced. 

It should also be noted that an alternative to choosing 
a model for the blurring aperture would be the estimation 
of the aperture directly from the data. However, in view 
of the extremely large quantization noise present in all the 
ERTS data examined and which will be discussed shortly, it 
seemed unlikely that this procedure would provide any better 
estimate of the blurring process. 

By re-examining the results of Table 4.1 and 4.2, it is 
apparent that the more data samples per blurring aperture 
width, the more effective the restoration procedure for a 
given set of constraints. For example in the extreme case 
where the blurring aperture is described by only one sample 
point, then regardless of the constraints on the restoration 
filter the resulting composite system point- spread function 


can contain no less than one sample point. Thus for this 
case it is clear that no resolution improvement is possible. 
However, if the blurring aperture could be sampled more often, 
increasing the number of sample points per blurring width, 
then a resolution improvement should be possible. 

This problem of having relatively few samples per 
blurring aperture width exits with the ERTS data. The 
scanner aperture has an effective ground diameter of approx- 
imately 80 meters. The data is collected so that a rectan- 
gular array of points is generated with an effective ground 
distance between adjacent vertical points in the array of 
79 meters and an effective ground distance between adjacent 
horizontal points in the array of 56 meters. Thus it is 
unlikely that there are more than 3 sample points within a 
circle of three a radius centered at the blurring aperture. 

In order to obtain significant resolution improvement, it is 
necessary to have a larger number of samples per blurring 
width than this. One way to obtain more samples is to 
interpolate additional values between the original data 
points. The idea is not to create information where no 
information originally existed buth rather to provide suf- 
ficient sample values to allow use of an effective restoration 
filter function. The procedure that will be used is to 
increase the number of points per blurring aperture until 
they correspond to one of the aperture functions for which 
the optimum restoration filter has been determined. The 



[ 


1 


I 


I 


I 


144 


The restoration can then be carried out and the final data 
set can be compressed back to the original size if desired. 

The use of interpolation allows the data set to be matched to 
the restoration function experimentally which is very desirable 
because precise values for the actual blurring function are 
not known. 

At the same time that the interpolation is performed to 
increase the number of points per blurring aperture, the 
scale factor differences of the vertical and horizontal image 
axes can be corrected. This scale factor difference arises 
from the fact that the effective ground distances between 
adjacent points along the vertical and horizontal axes are 
not equal. If these points are displayed as if these dis- 
tances were equal, then resulting distortion produces an 
image which is compressed vertically. Thus to correct this 
distortion, a somewhat greater degree of interpolation, or 
magnification, is required along the vertical axis. By 
magnifying the vertical axis, 79/56 with respect to the 
horizontal axis, this distortion should be corrected. 

A third-order polynominal interpolation function was 
used because of its relative smoothness compared to other 
forms, including trigonometric and sine functions. The 
restoration filter corresponding to a blurring aperture of 
radius of gyration of 5 is the largest that has been computed 
due to limitations imposed by round-off errors in the 
computations. In this blurring function, there are approx- 
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imately 45 samples (pixels) within -the 1% amplitude levels 
of the blurring function along each axis. Therefore, it is 
necessary to interpolate the ERTS data to obtain this same 
number of points in the estimated width of the ERTS aperture. 
As a first approximation, the width of the ERTS aperture was 
taken to be 3 as discussed above and accordingly the required 
interpolation factor is 45/3 - 15. In addition, to correct 
the geometric distortion, the ratio of the vertical to 
horizontal interpolation factors should be 79/56 ~ 1.411. 

The only integer ratio which satisfied both the interpolation 
factor for matching the effective data scanner aperture to a 
Gaussian blurring aperture with a radius of gyration of five 
and the geometric distortion correction interpolation ratio 
was 17/12 c 1.417. In other words, an interpolation of 17 
along the vertical axis and 12 along the horizontal axis 
would not only correct the geometric distortion but also 
approximately match the effective ERTS data scanner aperture 
to the restoration function computed for a Gaussian aperture 
with a radius of gyration of five. Several additional 
interpolation ratios were used to determine whether a better 
aperture match could be obtained. However, as will be 
explained later, this ratio provided the best performance. 

It appears likely that restoration functions computed 
for larger radii ox gyration then five would provide even 
better results when matched to the ERTS data by appropriate 
interpolation factors. However, because of the excessive 
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round-off errors encountered in the restoration function 
computation program, the largest radius of gyration which 
would provide reliable results was five. In addition, 
excessive interpolation of under sampled data may be desir- 
able "because of the effect of aliasing errors on the inter- 
polation accuracy. However, for the ratios chosen, this 
error did not appear significant. 
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5.3 Composite System Noise 

The composite system noise is defined as the total noise 
process including both the noise introduced by the ERTS data 
collection system and that introduced by the restoration pro- 
cess. The noise introduced by the ERTS data collection system 
arises from three primary sources: the sensors of the multi- 

spectral scanner, the signal conditioning electronics, and the 
A/D conversion process. The noise introduced by the resto- 
ration process arises principally from numerical round-off 
errors. 

The predominant data collection system noise source is 
introduced in the quantization stage of the A/D process. The 
data is quantized into six bits producing a signal dynamic 
range from 0 to 63. The total noise power, n 2 , introduced 
by uniform quantization may be approximated by the expression 


[3,6] 


n 


Av 

12 


(5-1) 


where Av represents the quantization increment, in this case 
it is 1 unit. Thus the total quantization noise power is 
approximately 

n 2 = J=r ~ 0.08333 volts 2 . (5-2) 


This noise process may be shown to have an essentially 
flat spectrum extending to many times the sampling r^queucy. 
A conservative estimate of the power spectrum density of 
the quantization noise can be obtained by assuming that it 
lies entirely within the band occupied by the signal. The 
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resulting one sided noise spectral density is given by 


$ (f) 
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f > 0 


(5-3) 


where f is the sampling frequency of the A/D process. The 

magnitude of f depends on the units of measurement being 

employed. For convenience it will be assumed that f =128 

s 

which is the number of samples used to define the restoration 
filters along each axis. With this convention the noise 
spectral density is then 

*nn lf) ^ 6TT28T = 1 ‘ 302 x 10 ‘ 3 ' 

Because the ERTS data does not occupy the full 64 level 
dynamic range defined by the quantization process, the 
quantization noise is made even more significant due to the 
reduced input signal-to-noise ratio. The ERTS data sets 
examined have a dynamic range of approximately 1 15 units, 
producing an input signal-to-noise ratio of less than 34 dB. 

In view of the fact that the restoration process will be 
enhancing the high frequency components of the ERTS data and 
from the assumption that the quantization noise spectral 
density is flat, the quantization noise presents a very 
serious problem to obtain significant restoration and a high 
signal-to-noise ratio in the deblurred image. 

The second source of composite system noise is introduced 
by computational round-off errors. These errors arise 
principally within the FORT Fast Fourier Transform algorithm. 

A double precision version of this algorithm is employed both 
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in the computation of the restoration function and in the 
application of this function to a blurred image. This 
algorithm is used to perform the spatial deconvolution by 
multiplication of the spectral components of the blurred 
image with the corresponding components of the restoration 
function and then talking the inverse Fourier transform of 
the product. By using this algorithm to compute the spectrum 
of one of the Gaussian blurring apertures previously de- 
scribed in Chapter 4 which should also be Gaussian shaped/ 

it was found that the magnitudes of the resulting spectrum 

-5 

never become less then approximately 10 when normalized 
with respect to the DC term. Consequently, this noise process 
is assumed to have a flat spectral density of approximately 

-5 

10 below the DC level of the data being corrected. Since 

the ERTS data sets examined typically had a DC value of 

approximately 30, the round-off error noise spectral density 

-4 2 

amplitude was assumed to be approximately 10 v /Hz. 

The severity of this problem is made evident by con- 
sidering the application of the classical "inverse" resto- 
ration filter to the image of Figure 4.106 as discussed in 
Section 4.4. Using the exact inverse one would expect to 
obtain a perfect restoration. However, because of the noise 
introduced by the round-off errors within the FORT algorithm, 
the actual restoration obtained consisted only of a random 
noise pattern with none of the original features being 
evident. 
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Since this noise source appears to be a function of the 
numerical accuracy of the computations within the FORT 
algorithm it is probable that it could be reduced to a 
negligible level by increasing the precision of the com- 
putations. However, it was not possible to employ more than 
double precision accuracy because of computer limitations. 

5.4 Resolution Enhancement of ERTS Data 

The ERTS data set to be enhanced had to satisfy two 
primary criteria. The first, and possibly the most important, 
was that the data set had to be free of any obvious noise so 
as to permit maximum restoration. The second was that the 
data set had to include a known topographical feature to 
provide a degree of "a priori" knowledge of the original 
image before blurring by the multispectral scanne. for com- 
parison to the restored results. Figure 5.1 and 5.2 show 
the ERTS data set selected. This data was taken from ERTS 
frame ID no. 1080-15192 collected on October 11, 1972, and 
reformatted as LARS run number 7^041900, Ch. 3, lines 916- 
1427, columns 986-1497. This area includes a portion of 
Washington, D.C. The topographical feature selected for 
enhancement was the Pentagon Building which appears at the 
center of Figure 5.1 and the upper center of Figure 5.2, a 
2X linear magnification of Figure 5.1. This feature was 
chosen not only because its shape is well defined but also 
because the luminance intensity of a single scan line passing 
through the center of the building for a scanner with 



Figure 5.1 ERTS Multispectral Scanner Data 
of Washington, D.C. Area with Pentagon at Center 
(LARS Run No. 72041900, Ch 3, Lines 916-1427, 
Columns 986-1497) 


Figure 5.2 2X Linear Enlargement of Figure 5.1 
Showing Pentagon in Upper Center 
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infinite resolution would appear as two closely spaced 
rectangular pulses with the level between the pulses re- 
presenting the "hole" in the center of the building. Thus 
by comparing the luminance signal for a single scan line 
passing through this area, a qualitative measure of the 
resolution enhancement is possible. 

Figure 5 . 3 shows the Pentagon and surroundings of Figure 
5.1 after the application of a third-order two-dimensional 
polynomial interpolation. As previously discussed, in order 
to correct the axis scale factor difference and also to 
approximately match the ERTS scanner aperture to the most 
powerful restoration functions computed for a Gaussian blur- 
ring aperture with a radius of gyration of five, the hori- 
zontal axis was interpolated or magnified by a factor of 12 
while the vertical axis was interpolated by a factor of 17. 
Figure 5.4 represents a 2X linear enlargement of the upper 
left quadrant of Figure 5.3. 

Figure 5.5 represents a plot of the luminance level for 
a scan line passing through the center of the Pentagon of 
Figure 5.3, i.e. data line 110. The left most side of the 
Pentagon is centered at the 60th sample, the "hole" at the 
85th sample, and the right most side at the 110th sample. 

The frequency spectrum corresponding to this line is shown 
in Figure 5.6. The information provided by these two graphs 
will form the primary basis for selecting the "best" resto- 
ration function of those computed in Chapter 4 for a Gaussian 



Figure 5.3 Polynomial Interpolation of 
Figure 5.1 Using an Interpolation Factor of 
17 for the Vertical Axis and 12 for the Horizontal 
Axis 
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Figure 5.4 2X Linear Enlargement of the Upper 
Left Quadrant of Figure 5.3 
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aperture with a radius of gyration of five. 

The restoration function selection procedure was as 
follows. Several restoration functions were computed based 
upon white noise spectral density amplitudes approximating 
the spectral density amplitude estimates of the system noise 
processes previously described. The data of Figure 5.3 was 
restored by each function. Then the graph of data line 110 
and the corresponding frequency spectrum were compared to 
those of Figure 5.5 and 5.6. The spectrum of line 110 for 
each restoration was examined to determine if the fundamental 
frequency components of Figure 5.5, approximately the first 
22 normalized frequency components of Figure 5.6, had been 
amplified without introducing higher frequency noise com- 
ponents of comparable magnitude. In addition, the graphs of 
data line 110 for each restoration were compared to that of 
Figure 5.5 to determine the degree of resolution enhancement 
by estimating the rise time associated with the luminance 
level change between the left most side and the "hole" of 
the Pentagon. 

Figure 5.7 and 5.8 represent the graphs of data line 

110 and its spectrum when the restoration function based on 

— 5 

flat noise spectrum $ nn (f) = 10 (Figure 4.67) is applied 
to the data of Figure 5.5. The spectrum shorn in Figure 5.8 
indicates that significant high frequency spectral components 
have been introduced during the restoration procedure. 

Since these high frequency components have the same general 
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Figure 5.7 Graph of Data Line 110 of Figure 5.5 
after Restoration with the Filter of Figure 4.67 



Figure 5.8 Frequency Spectrum of Data Line 110 
of Figure 5.7 


shape as the spectrum of the restoration filter, shown in 

Figure 4.65, it is reasonable to assume that these spectral 

components are the result of amplification of the system 

noise by the restoration filter. In other words, the assumed 

-5 

noise spectral density amplitude of 10 for the computation 
of this filter was much less than that of the composite 
system noise process. Figure 5.9 shows the result of this 
restoration to Figure 5.3. Figure 5.10 is a 2X linear 
magnification of the upper left quadrant of Figure 5.9. 

By increasing the assumed noise spectral density 

-4 

amplitude to 10 a new restoration functxon shown in Figure 
4.73 is obtained. Figure 5.11 and 5.12 show the graph of 
data line 110 and its corresponding spectrum when this 
restoration function is applied to Figure 5.3. It is clear 
that the high frequency noise components, shown in Figure 
5.12, are of less amplitude then those of Figure 5.8. It 
should again be noted that the shape of the noise spectrum 
in Figure 5.12 is determined by the spectrum of the resto- 
ration filter. Figure 5.13 represents the result of apply- 
ing the restoration function of Figure 4.73 to the data of 
Figure 5.3. Figure 5.14 is a 2X linear magnification of the 
upper left quadrant of Figure 5.13. 

By further increasing the noise spectral density 
-3 

amplitude to 10 , the restoration function of Figure 4.79 

is obtained. Figure 5.15 and 5.16 represent the graphs of 
data line 110 and its corresponding frequency spectrum when 


Figure 5.9 Restoration of Figure 5.3 with the 
Filter of Figure 4.67 



Figure 5.10 2X Linear Enlargement of the Upper 
Left Quadrant of Figure 5.9 




Figure 5.13 Restoration of Figure 5.3 with the 
Filter of Figure 4.73 
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the restoration function of Figure 4.79 is applied to the data 
of Figure 5.3. Although the high frequency noise components 
are still evident in Figure 5.16/ the effect on the data, as 
shown in Figure 5.15, appears negligible. Figure 5.17 re- 
presents the restoration of Figure 5.3 by the restoration 
function of Figure 4.79. Figure 5.18 is a 2X linear magni- 
fication of the upper left quadrant of Figure 5.17. 

Additional restoration functions were computed based on 
various other amplitudes and shapes of the noise Spectral 
density. These correspond to Figure 4.70, 4.73, 4.75, 4.83, 
4.87, and 4.91, and each was applied to the data of Figure 
5.3. However, no significant improvement was achieved com- 
pared to the results of Figure 5.17 and 5.18 obtained by the 
restoration function of Figure 4.79. It should also be noted 
that additional interpolation ratios were applied to the image 
of Figure 5.1 to determine if a better "match” could be 
obtained than that of 17 for the vertical axis and 12 for the 
horizontal axis. These ratios included 10/7, 14/10, and 
20/14. However, no improvement was obtained over the 17/12 
ratio for any of the restoration functions previously dis- 
cussed. 

A comparison of the graph of data line 110 for this 
restoration. Figure 5.15, with that of the original inter- 
polated image. Figure 5.5, shows that the rise time for the 
restoration as measured for the transition between the left 
most side of the Pentagon and the "hole” is approximately one 
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half that of the original. In addition, the sides of the 
Pentagon as shown in Figure 5,17 and 5.18 are displayed with 
fewer grey levels, the transition defining the sides is more 
sharp and the "hole" is somewhat larger than as shown in 
Figure 5.3 and 5.4, as would be expected by a resolution 
enhancement. There is also better delineation between the 
sides of the Pentagon and the topographical features imme- 
diately surrounding it in the restored images of Figure 5.17 
and 5.18. Also an examination of the other topographical 
features of Figure 5.3 and 5.4 shows that these features are 
displayed with greater contrast and are less amorphous in 

ft 

shape in the restored images of Figure 5.17 and 5.18. 

One of the most significant results of the restorations 
shown in Figure 5.17 and 5.18 is with regard to the data 
truncation error. For the classical "inverse" restoration 
filter, with no constraint of the radius of gyration of its 
point-spread function, this error would have propagated 
throughout the corrected data set substantially reducing 
image quality. As shown by the restorations in these figures, 
this error though noticeable at the extreme edges, becomes 
negligible after about 25 points into the data record from 
each edge. In addition to reducing the data truncation 
error, the constraint in the radius of gyration of the point- 
spread function of the restoration filter also produced a 
function which contained less than 30 points of significant \ 

j 

amplitude as compared to the classical "inverse" restoration I 

i 

j 

i 
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point-spread function which would have required as many points 
as the data record to be corrected, which for the corrections 
of Chapter 5 would have been 512 points. This reduction in 
the size of the restoration function should also increase 
data processing efficiency. 

5.4 Problems of ERTS Data Collection System 

As a consequence of this study, two significant problems 
were encountered that are inherent in the ERTS data collection 
system. The first problem is the poor srgnal-to-noise ratio 
of the ERTS data. This problem could be reduced by both 
reducing the quantization level increment and by increasing 
the dynamic range of the data prior to quantization. The 
second problem, the need for relatively large interpolation 
ratios to match the effective blurring aperture of the ERTS 
data system to the available restoration functions, could be 
minimized by increasing the data sampling rate. Since less 
interpolation would be required, both greater interpolation 
accuracy and increased processing efficiency could be 
achieved. By address ing both of these problems, more power- 
ful restoration functions may be applied to the data thus 
producing even greater resolution improvement. 


u ........ J J -1 1 1 - J 
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CHAPTER 6 

CONCLUSIONS AND SUGGESTIONS FOR FURTHER STUDY 
6.1 Conclusions 

This study has demonstrated that it is possible to design 
a restoration filter based upon minimizing the radius of 
gyration of the composite system point-spread function while 
constraining the radius of gyration of the point-spread 
function of the restoration filter and the noise power in the 
restored image. By assuming that the blurring aperture and 
composite system noise process are separable, a significant 
reduction in solution complexity was introduced. Transforma- 
tion of the defining spatial equations into the frequency 
domain led to a further simplification in that a system of 
linear differential equations results for second-order 
spatial weighting funcations in both the radii of gyration 
of the composite system and restoration filter point-spread 
function. Through application of well-known numerical 
solution techniques to this system of equations, a solution 
which satisfied the original constraints was obtained. 
Additional control over secondary oscillations in the com- 
posite system point-spread function was made possible by 
the introduction of an iterative technique. It was further 
shown that the resulting restoration filters satisfied the 
fundamental design criteria: reduction of the radius of 

gyration of the composite system point-spread function, and 
reduction of truncation error by constraining the radius of 
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gyration of the restoration filter point-spread function. 

One of the most significant results of this study was 
the successful application of this restoration technique to 
data which was blurred by an unknown function. By assuming a 
Gaussian model for the ERTS multispectral scanner aperture 
and employing a interpolation technique to match the effective 
scanner blurring aperture to the Gaussian aperture for which 
the most powerful restoration functions had been computed, a 
significant restoration of truncated ERTS data was obtained. 

In addition, because of the constraint on the radius of 
gyration of the restoration filter point- spread function, the 
truncation error was almost entirely limited to a small region 
near the edges of the image. Thus substantial restorations 
having a high signal-to-noise ratio and negligible truncation 
error were achieved. 

6.2 Suggestions for Further Study 

As a consequence of this investigation, several areas of 
additional study may be defined. First by compressing the 
enhanced interpolated data to its original scale, classifi- 
cation studies may be performed to determine the effect of 
resolution enhancement in classification accuracy. Secondly, 
by increasing the numerical precision of the restoration 
function computation program, more powerful restoration 
functions computed for blurring functions having more points 
per blurring aperture width than that for the Gaussian aper- 
ture with a radius of gyration of five could be applied to 


I 


1 


1 
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the ERTS data. Thirdly, the design of a restoration filter 
which would have a spatially variant constraint on the radius 
of gyration of the point-spread function of the restoration 
filter should provide greater resolution enhancement near the 
center of a truncated image while simultaneously controlling 
truncation error at the edges of the restoration than the 
spatially invariant constraint on the filter of this study. 
Fourthly, the development of a similar filter using finite 
difference equations rather than differential equations may 
reduce the numerical errors and resulting solution in- 
stabilities encountered with differential equation solution 
algorithms. Finally the analysis of blurring apertures for 
other than Gaussian functions should provide useful results, 
particularly for rectangular functions where linear motion 
blur correction is desired. 




I 
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APPENDIX A 

DERIVATION OF RELATIONSHIP BETWEEN DATA RECORD 
TRUNCATION ERROR AND RESTORATION FILTER POINT- 
SPREAD FUNCTION 

For clarity, the following analysis is based upon 
one-dimensional blurring and restoration functions. 

Figure A-l defines the parameters of the basic system under 
consideration. The blurring system output, y(t), can be 
expressed as 

CO 

y(t) = / x(T)h.(t-T)dT (A-l) 

— CO 

where x(t) is the original signal and h^tt) is the blurring 
function. A truncated record of y(t) could be defined as 

y T (t) - y(t) u(t-t Q ) 

® y(t) - y(t) u(t Q -t). (a-2) 

A 

Define z(t) to be the restored version of x(t) which is 

obtained from a correction of y^(t) as 

00 

z (t) = / h r (t-e)y T (g)d6 
— 00 


(A- 3) 


h b (t) h r (l) 


X(t) _ 

BLURRING 

y (t ) 


RESTORATION 

z(t)L 

7^ 

FUNCTION 


FUNCTION 

_ __l. — i. j 


Figure A-l: One Dimensional Composite Imaging 

System 
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which after substituting Eq. A- 3 into A- 2 . becomes 

00 

z(t) * / h r (t-e) [y(B)-y (3 )u(t Q -P)]d$ 

mOO 


00 

= / h r (t-6)y(0)d$ - e T (t) (A- 4 ) 

— 00 


for 


< t < » 


where (t) represents the error resulting from a 
truncation of y(t) and may be expressed as 


e (t) = / h (t-e) y(e)de 

1 —CO X 


(A- 5) 


for 

t. < t < ■ 

Define h rl (t-B), a duration limited form of h r (t-fS) so that 

h rl (t-6) = m(t-e) h r (t-B) (A-6) 

where m(t-B) is a monotomically decreasing weighting 
function such as expE-(t-B)*] so that 

m(0) - 1 (A-7) 


and 

liin (t-e)+“ m(t-3) - 0 


(A- 8) 


Prom Eq. AS and 6 define the truncation error resulting 
from a restoration with a duration limited form of h^. ( * ) as 
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0 

e Tl (t) » / h rl (t-0}y($)dB 

*— fca 


for 


t Q < t < ». 


(A-9 ) 


Prom Eg. A- 6 and 9, 


B {e T1 2 (t)} = E { / m(t-e)h v (t-(3)y((3)aB 


/ m(t-a)h (t-a)y (a)da> 

— oo 


= E { // m(t-$)m(t~a)h r (t-B)h r (t-a)y(6)y (a)dadp} 

“CO 

(A- 10) 

which because of the shape of the weighting function. 
Eg. A- 7 and 8, may be upper bounded 


'0 


E {e Tl 2 (t)} <_ m 2 (t)E{ // h^ {t-8 Jh^ (t-a )y (8 )y (a)dadg) . 


(A- 11) 


Assuming that y ( • ) is a sample function from a wide-sense 

stationary random process, then Eq. A-ll may be rewritten 

as . 

t 0 

E {e-, 2 ^} <_ m 2 (t)// h„(t-P)h (t-a)R . (a-B)dad0 

* — oo J; J 


< m 2 (t) R. ( 0 ) [ / h (t-&)d&] 2 . (A-12) 

¥ —CO 
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Since there is a constraint on the energy of the 
restoration function, the right side of Eq. A- 12 is finite 
for all t. Taking the limit with respect to t of both 
sides of Eq. A-12, 

lim t ^ E{e T1 2 (t)} - 0, (A-13) 


since 


lim^^ m 2 (t) - 0, from Eq. A- 8 . 


Thus as the distance from the truncation point approaches 
infinity, the variance of the truncation error produced by 
a duration limited restoration function approaches zero. 
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APPENDIX B 

DERIVATION OF EQUATION 2-21 


The gradient of Eq. 2-12 

CO 

VI = / [A X (U,S) + . 

— 00 

CO 

+ / [A 2 (u,z) 

— 00 

CO 

+ X 2 / [A 3 (u,z) 

— CO “ 

CO 

+ X 3 / [A 4 (u,i) 

— oo 


may be written as 
L£(u,z}]h r (z)dz 

+ A 2 (u / z}]h r (z)dz 

+ A 3 (u r z)]h r (z)dz 

+ A^ (u, z) ]h r (z)dz . (B-l) 


Substituting Eq. 2-14, 2-17 through 2-20 into Eq. B-l, 

CO 

/ [A^ (u, z ) + X 2 A 3 (u,z) + A 3 A 4 {u,z)]h r (z)dz 

—CO 

CO 

+ X. / s (h)8 (u-z )h (z )dz = 0 
X —00 


or 


—00 


/ [A x (u,£) + X 2 A 3 (u,z) + X 3 A 4 (u,z)]h r (z)dz 


+ X 1 s(u)h r (u) « 0, 


(B-2 ) 
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APPENDIX C 

DERIVATION OF EQUATIONS 2-25, 27, 29, and 31 
Substituting Eq. 2-13 into 2-24, 

CO 

= /// w{v)h b {v-z)h b Cv-u)dv 


i -j2irfu e j2irvz d ~ d ~ 


(c-l) 


Introducing a change of variable into Eq. C-l where 


v - u = a 


(C-2) 


B 1 (f,v) = /// w(v)h b (v-z')h b (a)e 

mCO 

j 2 ttvz - _- 

■ dv da dz 


-j27rf (v-a) 


e- 


= H^(£) //w(v)h b (v-s)e e^ 27r ^^dv dz 

(C-3) 


and introducing another change of variable into Eq. C-3 
where 

V - z = 3 ( C— 4 ) 


B x (f,v) = H*(I) // w{v)h b (3)e 


« '^~j 2ir f v ej2 7TV ^V“B ) d ^ d 3 


= ^(fJ^Cv) / w (v ) e 


-^-j27TV(f-V) d - 


* _ 


- ^ (f )H^(v)W(f-v), 


(C-5 ) 
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Following a similar development, 

B 2 (f,v) = \S B 2 (5,5)e-3 2,r55 a :i2 ” 5 u5ai (C-6) 

—00 

and substituting Eg. 2-14 into Eg. C-6, 

B 2 (f,v) = // s(u)6 (u-z)e”^ 2ir ^e^ 27r ^^dudz 

—CO 

= / S(5)e' j2,li(? ‘' :) d5 


= S(f-v). (C-7) 

B-(f,V) = // A 3 (u f z)e" j2,,r5u e j27TVZ dudz (C-8) 

—CO 

which after substituting Eg. 2-15 into Eg. C-8 becomes, 

B 3 (f,v) = // (z-u)e"'^ 2lT ^ u e^ 2ir ^ Z dudz. (C-9) 

— 00 

Introducing the change of variable, 

z - u = a (C-10) 

into Eg . C-9 , 


b <f,S) = 

— DO 

= »* W) 7 e- j2lr5<f -^>aS 

nn -L 

= (v)6(f-v). 

nn 


(C-ll) 
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And finally, 

B 4 (f,v) = j"/ A 4 (U,S) e -3 2 * f V 2 ” 5 fl5az 


(C-12) 


Substituting Eg. 2-16 into Eg. C-12, 


V f ' v) = — r~ 2 - 

4 TT 2 f „ 2 
H 


// h^'' (v-zjhj^ (v-u)dv 


= e 


- j 2»IU a j2»v5 d5d ; > 


(C-13 ) 


By introducing a change of variable into Eg. C-13 where 


v - z = a 


(C-14) 


oo co 


B 


a (f,5, = // { I h- (5)e^ 2 ”“a5 

4 IT f„ -*> - 00 

£i 


(C-15) 


and recognizing that 

/ hT' (H)e" j2irVa da = / {h b {a)*6~{5)}e” j ™da , (C-16) 

—CO D —CO 

and using the convolution property of the Fourier 
transform, 

/ h^'CoOe ^ 2lT ^da = - 4'!r 2 v 2 ^(v) . 


(C-17 ) 


Substituting Eg. C-17 into C-15, 


B 4 (f,v) = 

H 


( v ) //h b (v-u)e"^ 27r ^ U e^ 2lI ^ V dudv (c-18) 


and introducing another change of variable where 
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v - u = P , {C-19 ) 

Eq. C-18 becomes, 

b.(5,v) = -ll // 
f H 2 -» 

= jL |h. <v)| 2 ; e - 32 ^tf-v ) d5 

f H 
_ 2 

= |B b (v).f 2 6(f-v). (C-20) 

f TT 
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APPENDIX D 

DERIVATION OF EQUATION 2-92 
Substituting Eg. 2-47 and 50 into 2-91, 


* . • • • i 

O. r,i* * i 


71 = / [H (f )H (v ){ — -V ) 

X _4 bx x' DX x' ^ 2 X x' 

w 


+ c 6 {f -v ) - 


X 


4irx 


6"(v ) + cl(v-i)) 


w 


+ X, {- 


^ 4H 2 * 2 


6"" (f -v ) + d 8(f -v ) - 

' X X X X X 


1 6"(v -f ) 


4ir 2 x 2 x x 


+ d x 5 (VV }lH rx^x +2 ^nnx ( V 


or 


+ Vrr* WI'W ” °' 


(D-l } 


VI = 
x 


Be 

4 t£ *> _di 

2 *V df x 2 


[i k Cf x )H rx (f x )] 


+ 2c |h (f )| 2 h (f ) - 
x 1 bx x 1 rx x 


lx dr _ * 

— H (f j 

2-- 2 « 2 rx v x ; 


2tt 2 x df 

S X 


f x 2 


+ 2d x *W f x> + 2 + x 3x 




(D-2) 
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or 


VI- - - 'VSsl + 2^ x (f x )i^W x ) 


X 


2tt 2 x 2 

w 


.K 


+ VWV 1 - i^T OT 


Ik 


+ wy°' 


(D-3) 


or 


VI « 
x 


_ i ! ^ ( V l + — ]iC (0 


2 tt 2 x 2 
w 


2 u 2 x s 2 131 X 


*kx (X x^x (£ x^ , f % foV’kV 

A.* ** x 2 tt 2 x 2 


* *W 


w 


- 2 (-To +X (■ } H (£ ) 2 + A„ $ (f ) + d )] 

^ lc x Vf„ J J1 VV' 2x nnx v x' x 

Ek 


H (f ) “ Of 
rx' x 


(D-4 ) 
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or 


VI 


x 


[ 1 ( f X l |2y -S ‘V.’ (£ 


2ir 2 x 2 x 2 
w s 


rx x 


W f * )H bx (f *V , f , 

Z~2 


H X 


w 


1 * 
^ r tt 


- 4 ^ a V {IC x W 3x ( V f Hx> 21 


2* 2x w 2 


lyyl' + ^WV + VJW 


- 0 , 


(D-5) 


which may be written in normalized form as 


VI x = ^ (f x> + 2 ^x' f x )I £x tf x ,x s 2 ^x lf x ) 

/ nwviv - *ix x w 2] + vii'VK'v 

- 4* : V { f c x W 3x ' f x/ f HX 


+ *2x*nnx (f x> + a x }!H rx tf x>' /[ i H bx (f x , l 2x S 


2 
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APPENDIX E 

DERIVATION OP EQUATION 2-93 


From Eq. 2-37, 2-50, and 2-75, 

CO 


K 


lx 


< S x-V“™'V H rx<V a V f X 


+ d // 6 (f -v ) H* (f ) H (V )dv df 
x v x x' rx x x' rx v x' x x 


4tt 2 x_ 2 -L H rx t£ x ,H rx (£ x )af x 


s 


+ a x i l H rx (£ x>l 2a£ x- 


(E-l) 


From Eq. 2-94, 


E"(f ) 
rx x 


E" (£ ) + j H'r (f ) 
rrx x J rix x 


(E-2) 


and 


H <f ) 
rx x 


H (f ) - j H . (f 1 
rrx x' J rrx x* 


(E— 3 ) 


Multiplying Eq. E-2 and E-3, 


H'' " (f ) H (f ) = E" (f ) H (f } 

rx x rx x rrx x rrx x 


+ H"r (f )h . (f ) + j [H"r (f ) h (f ) 

rxx x rix x J rrx x rrx x 


- E" (f ) H . (f )]. 
rrx x rix x * 


(E-4 ) 
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Substituting Eq. E-4 into E-l, 


K 


lx 


i-t J(H'' (f )H Cf ) 


47r 2 x 2 rrx x' rrx x 

s 


+ H rix <f X > H ri X < f x > ,af X + 3_/ < H ri x < f x > 


H {f )-ff"* (f j H . (f )}df ] + d / J H (f )[ z df. 
rrx' x 7 rrx x rix x' x J x L 1 rx' x' 1 x 


Since h (x) is assumed to be a real function/ then 
r XI 


(E-5) 


and 


H"£„(f ) is an even function 

3T3TX X 


E TT (f ) is an odd function, 
rrx x 


Thus the second integral in Eq. E-5 is zero/ and 


K 


lx 


f EH'"* (f ) H {f )+H*? (f )H . ( f ) ] df 

J u rry v v * rrv ri v x rix X 


. 2 2 J rrx' x' rrx x' rrx' x' rix x 

4tt x - 00 


X 


+ d / I h (f )j 2 df , 
x J 1 rx x 1 v 

CO 


(E-6) 
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APPENDIX P 

DERIVATION OF EQUATION 2-119 


Substituting 


rx 


"1 0=1 

(T) = ? {H rx (f x )} = / H rx Cf x )e:]2lTfxTdf x (F " 1) 


into Eq. 2-116 r . 

I m = “2 // w v (x)m^(x)h hv (x-T) / H_ v (f v )e^ irf x T df v dTdx . 


rx x 


x 


Introducing the change of variables into Eq. F-2, 
x - T = a , 


(P-2) 

CF-3) 


then 

T m m " 2 / \x |f x > // w x U)m x U)h bx (a)e3 2,rf x (x - r!) 

«D0 — OO 

“ " 2 l H rx !f x )H b x (f x ) 7 w x te)m x (x)e :]2lrf x x axdf 

By making the substitution, 

w x (x) = ? tw x (v x )} = 7 W x (v x )e^ 27rv x x dv x 

— CO 


dxdadf 

A 


. (F— 4 ) 

X 


(F-5) 
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into Eq. 

F-4 , 

50 


I 

-2 / 

H 

m 

—CO 

rx 


dxdv 

df 


x x 


" ~ 2 i H rx (f * )B b x (f x >_/ VVVVV^-V tF “ 6) 

Finally Eq. F-6 may be written in the form of Eq. 2-118, 

CO 

I_ = / B_ (f )H (f )df (F-7 ) 

m 1 5x x rx x x 


where 


E 5x (f x> - - 2 ^x (f x> [ W^'W^x 


(F-8) 
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APPENDIX G 

TABULATION OP RESTORATION POINT- SPREAD FUNCTIONS 
OF FIGURES 4.67, 4.70, 4.73, 4.76, 4.79, 
4.83, 4.87, and 4.91 
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Figure 4.79 
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